home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / public / SciAn / src / ScianHwuFiles.c < prev    next >
C/C++ Source or Header  |  1994-08-01  |  153KB  |  5,299 lines
  1. /*******************************************************************************/
  2. /********************************STF READER ************************************/
  3. /*******************************************************************************/
  4. /*
  5.   10/26/1992
  6.   ScianHwuFiles.c
  7.   ReadSTFFile routine
  8.  
  9.   Tzong-Yow Hwu
  10.  
  11.   The file format is:
  12.   
  13.       RANK <rank>
  14.       DIMENSIONS <xdim> <ydim> ....
  15.       SCALAR|VECTOR <nComponents>
  16.       [ BOUNDS     <xmin> <xmax> <ymin> <ymax> .... ] {required when no grid}
  17.       [ INTERLACED|NONINTERLACED ]  { default is INTERLACED }  
  18.       [ ORDER COLUMN|ROW ]          { defalut is COLUMN major }
  19.       [ NAME <datasetName> ]
  20.       [ TIME <time> ]
  21.       DATA
  22.       <values> ......
  23.       .
  24.       .
  25.       .
  26.       VECTOR <nComponents>
  27.       [INTERLACED|NONINTERLACED] { takes previous spec as default }
  28.       [ ORDER COLUMN|ROW ]       { takes previous spec as default }
  29.       GRID
  30.       <values> ......
  31.       .
  32.       .
  33.       .
  34.       END    { the default is reset to INTERLACED and COLUMN major when 
  35.            anther datasets begin }
  36.       .
  37.       .
  38.       .
  39.       .
  40.  
  41.   Grid is optional for any dataset, if a grid is not present then regular
  42.   dataform is assumed.  The order of grid and data for a dataset is 
  43.   arbitrary in the data file.  However, rank and dimensions are required 
  44.   before either value is processed and the rank must precede dimensions. 
  45.   If no grid is specified, then bounds is required.  
  46.   If no order is specified for a set of values, then the default of 
  47.   column major is assumed.
  48.   If no dataset name is specified for the dataset, then it inherits from the
  49.   previous dataset in the file.  If this happens to the first dataset, then
  50.   it takes the name of the file.
  51. */
  52. #include "Scian.h"
  53. #include "ScianTypes.h"
  54. #include "ScianArrays.h"
  55. #include "ScianIcons.h"
  56. #include "ScianWindows.h"
  57. #include "ScianObjWindows.h"
  58. #include "ScianVisWindows.h"
  59. #include "ScianVisObjects.h"
  60. #include "ScianControls.h"
  61. #include "ScianColors.h"
  62. #include "ScianDialogs.h"
  63. #include "ScianFiles.h"
  64. #include "ScianFileSystem.h"
  65. #include "ScianLists.h"
  66. #include "ScianPictures.h"
  67. #include "ScianErrors.h"
  68. #include "ScianTimers.h"
  69. #include "ScianDatasets.h"
  70. #include "ScianFilters.h"
  71. #include "ScianTextBoxes.h"
  72. #include "ScianTitleBoxes.h"
  73. #include "ScianButtons.h"
  74. #include "ScianSliders.h"
  75. #include "ScianScripts.h"
  76. #include "ScianIDs.h"
  77. #include "ScianVisContours.h"
  78. #include "ScianStyle.h"
  79. #include "ScianSpaces.h"
  80. #include "ScianMethods.h"
  81. #include "ScianObjFunctions.h"
  82. #include "ScianHwuFiles.h"
  83.  
  84. #define STF_MAXHEADER 9
  85. #define STF_MAXNAMELEN 256
  86.  
  87. /* field type */
  88. #define STF_SCALAR 0
  89. #define STF_VECTOR 1
  90.  
  91. /* vector field data permutation type */
  92. #define STF_NONINTERLACED 0
  93. #define STF_INTERLACED 1
  94.  
  95. #define STF_COLUMNMAJOR 0
  96. #define STF_ROWMAJOR 1
  97.  
  98. /* for datasetType */
  99. #define STF_ERROR -1
  100. #define STF_END  0
  101. #define STF_DATA 1
  102. #define STF_GRID 2
  103.  
  104. /* for specFlag */
  105. #define STF_ISTYPE 1
  106. #define STF_ISDIMS 2
  107. #define STF_ISBOUNDS 4
  108. #define STF_ISNAME 8
  109. #define STF_ISORDER 16
  110. #define STF_ISRANK 32 
  111. #define STF_ISPERMUTATION 64 
  112. #define STF_ISTIME 128
  113.  
  114. typedef struct  {
  115.      int fieldType;   /* SCALAR or VECTOR */
  116.      long *dimensions; /* xdim, ydim, zdim .... */
  117.      real *bounds;  /* bounds array */
  118.      int nBounds; /* number of bounds provided by the user if any */
  119.      char datasetName[STF_MAXNAMELEN];
  120.      int order; /* COLUMNMAJOR(default) or ROWMAJOR */
  121.      int rank;
  122.      int permutation; /* INTERLACED or NONINTERLACED */
  123.      int nComponents;  /* number of components for vector data set */
  124.      real time;
  125.      } STF_spec;
  126.  
  127. #ifdef PROTO
  128. static int parseHeader(FILE *inFile, int *specFlag, STF_spec *datasetSpec)
  129. #else
  130. static int parseHeader(inFile, specFlag, datasetSpec)
  131. FILE *inFile;
  132. int *specFlag;
  133. STF_spec *datasetSpec;
  134. #endif
  135. {
  136.   int datasetType = STF_END; /* DATA, GRID, ERROR, or END */
  137.   int headerLine = 0;
  138.   char dirStr[360]; 
  139.   int counter;
  140.   int c;
  141.  
  142.   for (SkipBlanks(inFile); (c = getc(inFile)) >= 0 && headerLine < STF_MAXHEADER;
  143.        SkipBlanks(inFile))
  144.   {
  145.     ungetc(c, inFile);
  146.  
  147.     if (1 == fscanf(inFile, "%s", dirStr))
  148.       {
  149.         if (dirStr[0] == '#')
  150.           {
  151.             /* a comment */
  152.             ReadLn(inFile);
  153.             continue;
  154.           }
  155.  
  156.         if (0 == strcmp2(dirStr, "RANK"))
  157.           {
  158.             if ( *specFlag & STF_ISRANK )
  159.               {
  160.                 FileFormatError("ReadSTFFile", "Error reading header: double RANK specifications.");
  161.                 datasetType = STF_ERROR;
  162.                 break;
  163.               }
  164.             SkipBlanks(inFile);
  165.             if ( 1 != fscanf( inFile, "%d", &(datasetSpec->rank) ) )
  166.               {
  167.                 FileFormatError("ReadSTFFile", "Error reading header: RANK value format error.");
  168.                 datasetType = STF_ERROR;
  169.                 break;
  170.               }
  171.             if ( datasetSpec->rank < 1  )
  172.               /* rank must be greater than or equal to 1 */
  173.               {
  174.                 FileFormatError("ReadSTFFile", "Error reading header: invalid RANK value.");
  175.                 datasetType = STF_ERROR;
  176.                 break;
  177.               }
  178.              ReadLn(inFile);
  179.              headerLine++;
  180.              *specFlag |= STF_ISRANK;
  181.              continue; 
  182.            }
  183.  
  184.         if (0 == strcmp2(dirStr, "TIME"))
  185.           {
  186.             if ( *specFlag & STF_ISTIME )
  187.               {
  188.                 FileFormatError("ReadSTFFile", "Error reading header: double TIME specifications.");
  189.                 datasetType = STF_ERROR;
  190.                 break;
  191.               }
  192.             SkipBlanks(inFile);
  193.             if ( 1 != fscanf( inFile, "%g", (float *) &(datasetSpec->time) ) )
  194.               {
  195.                 FileFormatError("ReadSTFFile", "Error reading header: invalid TIME value.");
  196.                 datasetType = STF_ERROR;
  197.                 break;
  198.               }
  199.  
  200.              ReadLn(inFile);
  201.              headerLine++;
  202.              *specFlag |= STF_ISTIME;
  203.              continue; 
  204.            }
  205.  
  206.          if (0 == strcmp2(dirStr, "SCALAR"))
  207.            {
  208.              if ( *specFlag & STF_ISTYPE )
  209.                {
  210.                  FileFormatError("ReadSTFFile", "Error reading header: double specifications of SCALAR or VECTOR header line.");
  211.                  datasetType = STF_ERROR;
  212.                  break;
  213.                }
  214.              headerLine++;
  215.              ReadLn(inFile);
  216.              datasetSpec->fieldType = STF_SCALAR;
  217.              datasetSpec->nComponents = 0;
  218.              *specFlag |= STF_ISTYPE;     
  219.              continue;
  220.            }
  221.  
  222.          if (0 == strcmp2(dirStr, "VECTOR"))
  223.            {
  224.              if ( *specFlag & STF_ISTYPE )
  225.                {
  226.                  FileFormatError("ReadSTFFile", "Error reading header: double specifications of SCALAR or VECTOR header line.");
  227.                  datasetType = STF_ERROR;
  228.                  break;
  229.                }
  230.  
  231.              SkipBlanks(inFile);
  232.              if ( 1 != fscanf(inFile, "%d", &(datasetSpec->nComponents)) )
  233.                {
  234.                  FileFormatError("ReadSTFFile", "Error reading header: VECTOR component number format error.");
  235.                  datasetType = STF_ERROR;
  236.                  break;
  237.                }
  238.              if ( datasetSpec->nComponents < 1 )
  239.                {
  240.                  FileFormatError("ReadSTFFile", "Error reading header: invalid VECTOR component value.");
  241.                  datasetType = STF_ERROR;
  242.                  break;
  243.                }
  244.  
  245.              headerLine++;
  246.  
  247.              ReadLn(inFile);
  248.              datasetSpec->fieldType = STF_VECTOR;
  249.              *specFlag |= STF_ISTYPE;     
  250.              continue;
  251.            }
  252.  
  253.          if (0 == strcmp2(dirStr, "INTERLACED"))
  254.            {
  255.              if ( *specFlag & STF_ISPERMUTATION )
  256.                {
  257.                  FileFormatError("ReadSTFFile", "Error reading header: double INTERLACED/NONINTERLACED specifications.");
  258.                  datasetType = STF_ERROR;
  259.                  break;
  260.                }
  261.  
  262.              headerLine++;
  263.              ReadLn(inFile);
  264.              datasetSpec->permutation = STF_INTERLACED;
  265.              *specFlag |= STF_ISPERMUTATION;
  266.              continue;
  267.            }
  268.  
  269.          if (0 == strcmp2(dirStr, "NONINTERLACED"))
  270.            {
  271.              if ( *specFlag & STF_ISPERMUTATION )
  272.                {
  273.                  FileFormatError("ReadSTFFile", "Error reading header: double INTERLACED/NONINTERLACED specifications.");
  274.                  datasetType = STF_ERROR;
  275.                  break;
  276.                }
  277.  
  278.               headerLine++;
  279.               ReadLn(inFile);
  280.               datasetSpec->permutation = STF_NONINTERLACED;
  281.               *specFlag |= STF_ISPERMUTATION;
  282.               continue;
  283.            }
  284.  
  285.          if (0 == strcmp2(dirStr, "NAME"))
  286.            {
  287.              if ( *specFlag & STF_ISNAME )
  288.                {
  289.                  FileFormatError("ReadSTFFile", "Error reading header: double NAME specifications.");
  290.                  datasetType = STF_ERROR;
  291.                  break;
  292.                }
  293.  
  294.              headerLine++;
  295.              SkipBlanks(inFile);
  296.  
  297.              if ( 1 != fscanf(inFile, "%s", datasetSpec->datasetName))
  298.                {
  299.                  FileFormatError("ReadSTFFile", "Error reading header: NAME specificaton format error.");
  300.                  datasetType = STF_ERROR;
  301.                  break;
  302.                }
  303.              ReadLn(inFile);
  304.              *specFlag |= STF_ISNAME;
  305.              continue;
  306.            }
  307.  
  308.          if (0 == strcmp2(dirStr, "ORDER"))
  309.            {
  310.              if ( *specFlag & STF_ISORDER )
  311.                {
  312.                  FileFormatError("ReadSTFFile", "Error reading header: double ORDER specifications.");
  313.                  datasetType = STF_ERROR;
  314.                  break;
  315.                }
  316.  
  317.              headerLine++;
  318.              SkipBlanks(inFile);
  319.  
  320.              if ( 1 != fscanf(inFile, "%s", dirStr))
  321.                {
  322.                  FileFormatError("ReadSTFFile", "Error reading header: ORDER specification format error.");
  323.                  datasetType = STF_ERROR;
  324.                  break;
  325.                }
  326.              if ( 0 == strcmp2(dirStr, "COLUMN") )
  327.                {
  328.                  datasetSpec->order = STF_COLUMNMAJOR;
  329.                }
  330.              else if ( 0 == strcmp2(dirStr, "ROW") )
  331.                {
  332.                  datasetSpec->order = STF_ROWMAJOR;
  333.                }
  334.              else
  335.                {
  336.                  FileFormatError("ReadSTFFile", "Error reading header: invalid ORDER specification.  It should be COLUMN nor ROW.");
  337.                  datasetType = STF_ERROR;
  338.                  break;
  339.                }
  340.  
  341.              ReadLn(inFile);
  342.              *specFlag |= STF_ISORDER;
  343.              continue;
  344.            }
  345.  
  346.          if (0 == strcmp2(dirStr, "DIMENSIONS"))
  347.            {
  348.              if ( !(*specFlag & STF_ISRANK) )
  349.                {
  350.                  FileFormatError("ReadSTFFile", "Error reading header: RANK missing befor DIMENSIONS");
  351.                  datasetType = STF_ERROR;
  352.                  break;
  353.                }
  354.                  
  355.              if ( *specFlag & STF_ISDIMS )
  356.                {
  357.                  FileFormatError("ReadSTFFile", "Error reading header: double DIMENSIONS specifications.");
  358.                  datasetType = STF_ERROR;
  359.                  break;
  360.                }
  361.  
  362.              headerLine++;
  363.  
  364.             /* allocate storage for dimensions*/
  365.             if ( !(datasetSpec->dimensions = (long *) Alloc(datasetSpec->rank * sizeof(long))) )
  366.               {
  367.                 FileFormatError("ReadSTFFile", "Unable to allocate memory");
  368.                 datasetType = STF_ERROR;
  369.                 break;
  370.               }
  371.  
  372.              *specFlag |= STF_ISDIMS;
  373.  
  374.              /* read all the dimensions */
  375.              for(counter = 0; counter < datasetSpec->rank; counter++)
  376.                 {
  377.                   SkipBlanks(inFile);
  378.                   if ( 1 != fscanf(inFile, "%ld", datasetSpec->dimensions + counter) )
  379.                     {
  380.                       FileFormatError("ReadSTFFile", "Error reading header: DIMENSIONS value format error. ");
  381.                       datasetType = STF_ERROR;
  382.                       break;
  383.                     }
  384.                 }
  385.  
  386.              if (datasetType == STF_ERROR)
  387.                break;
  388.  
  389.              ReadLn(inFile);
  390.              continue;
  391.            }
  392.  
  393.          if (0 == strcmp2(dirStr, "BOUNDS"))
  394.            {
  395.              if ( *specFlag & STF_ISBOUNDS )
  396.                {
  397.                  FileFormatError("ReadSTFFile", "Error reading header: double BOUNDS specifications.");
  398.                  datasetType = STF_ERROR;
  399.                  break;
  400.                }
  401.  
  402.              headerLine++;
  403.  
  404.              /* allocate storage for bounds*/
  405.              if ( !(datasetSpec->bounds = (real *) Alloc(sizeof(real))) )
  406.                {
  407.                  FileFormatError("ReadSTFFile", "Unable to allocate memory.");
  408.                  datasetType = STF_ERROR;
  409.                  break;
  410.                }
  411.  
  412.              *specFlag |= STF_ISBOUNDS;
  413.  
  414.              /* read all the bounds */
  415.              for( counter = 0, SkipBlanks(inFile);
  416.                   (c = getc(inFile))!= '\n' && c >= 0;
  417.                   counter++, SkipBlanks(inFile)
  418.                 )
  419.                 {
  420.                   ungetc(c, inFile);
  421.                   if ( !( datasetSpec->bounds = (real *) Realloc(datasetSpec->bounds, (counter+1)*sizeof(real))) )
  422.                     {
  423.                       FileFormatError("ReadSTFFile", "Unable to allocate memory.");
  424.                       datasetType = STF_ERROR;
  425.                       break;
  426.                     }
  427.                   if ( 1 != fscanf(inFile, "%f", (float *) (datasetSpec->bounds + counter)) )
  428.                     {
  429.                       FileFormatError("ReadSTFFile", "Error reading header: BOUNDS value format error.");
  430.                       datasetType = STF_ERROR;
  431.                       break;
  432.                     }
  433.                  }
  434.              
  435.              if (datasetType == STF_ERROR)
  436.                /* this break gets out of the outer for loop */
  437.                break;
  438.  
  439.              /* the number of bounds must be an even number */
  440.              if ( counter == 0 || counter % 2 != 0 )
  441.                {
  442.                  FileFormatError("ReadSTFFile", "Error reading header: BOUNDS value format error.");
  443.                  datasetType = STF_ERROR;
  444.                  break;
  445.                }
  446.  
  447.              datasetSpec->nBounds = counter;
  448.              continue;
  449.            }
  450.  
  451.          if (0 == strcmp2(dirStr, "DATA"))
  452.            {
  453.              headerLine++;
  454.              ReadLn(inFile);
  455.              datasetType = STF_DATA;
  456.              break;
  457.            }
  458.  
  459.          if (0 == strcmp2(dirStr, "GRID"))
  460.            {
  461.              headerLine++;
  462.              ReadLn(inFile);
  463.              datasetType = STF_GRID;
  464.              break;
  465.            }
  466.  
  467.          if (0 == strcmp2(dirStr, "END"))
  468.            {
  469.              ReadLn(inFile);
  470.              break;
  471.            }
  472.  
  473.          FileFormatError("ReadSTFFile, unknown header line specification", dirStr);
  474.          datasetType = STF_ERROR;
  475.          break;
  476.        }
  477.      else
  478.        {
  479.          /* a blank line*/
  480.          ReadLn(inFile);
  481.          continue;
  482.        }
  483.   }
  484.  
  485.   if ( datasetType == STF_END )
  486.     {
  487.       if ( headerLine != 0 )
  488.         {
  489.           FileFormatError("ReadSTFFile",
  490.           "Error reading header: Incomplete header specification for dataset.");
  491.           datasetType = STF_ERROR;
  492.         }
  493.       else
  494.         {
  495.           return datasetType; 
  496.         }
  497.     }
  498.  
  499.   if ( 
  500.        (datasetType == STF_ERROR) || 
  501.  
  502.        /* rank and dimensions are always required */
  503.        ( !(*specFlag & STF_ISRANK) || !(*specFlag & STF_ISDIMS) ) ||
  504.  
  505.        /* for grid data set, fieldtype is required */
  506.        /* and its fieldtype must be vector */
  507.        ( (datasetType == STF_GRID) && (!(*specFlag & STF_ISTYPE) || (datasetSpec->fieldType != STF_VECTOR) ) ) ||
  508.        
  509.        /* for field data set, fieldtype is required */
  510.        /* and for vector dataset, permutation is required */
  511.        ( (datasetType == STF_DATA) && (!(*specFlag & STF_ISTYPE)) ) 
  512.      )
  513.     {
  514.       switch ( datasetType )
  515.         {
  516.           case STF_ERROR:
  517.             /* error message was broadcasted */
  518.             break;
  519.           case STF_DATA:
  520.             if ( !(*specFlag & STF_ISRANK) )
  521.               {
  522.                 FileFormatError("ReadSTFFile",
  523.                 "Error reading header: RANK missing before DATA.");
  524.               }
  525.             else if ( !(*specFlag & STF_ISDIMS) )
  526.               {
  527.                 FileFormatError("ReadSTFFile",
  528.                 "Error reading header: DIMENSIONS missing before DATA.");
  529.               }
  530.             else if ( !(*specFlag & STF_ISTYPE) )
  531.               {
  532.                 FileFormatError("ReadSTFFile",
  533.                 "Error reading header: SCALAR/VECTOR missing before DATA.");
  534.               }
  535.             break;
  536.           case STF_GRID:
  537.             if ( !(*specFlag & STF_ISRANK) )
  538.               {
  539.                 FileFormatError("ReadSTFFile",
  540.                 "Error reading header: RANK missing before GRID.");
  541.               }
  542.             else if ( !(*specFlag & STF_ISDIMS) )
  543.               {
  544.                 FileFormatError("ReadSTFFile",
  545.                 "Error reading header: DIMENSIONS missing before GRID.");
  546.               }
  547.             else if ( !(*specFlag & STF_ISTYPE) )
  548.               {
  549.                 FileFormatError("ReadSTFFile",
  550.                 "Error reading header: VECTOR missing before GRID.");
  551.               }
  552.             else if ( datasetSpec->fieldType != STF_VECTOR )
  553.               {
  554.                 FileFormatError("ReadSTFFile",
  555.                 "Error reading header: GRID must be vector type");
  556.               }
  557.             break;
  558.           default:
  559.               FileFormatError("ReadSTFFile",
  560.               "Error reading header: header format error");
  561.             break;
  562.         }
  563.       datasetType = STF_ERROR;
  564.     }
  565.   return datasetType;
  566. } /* parseHeader */
  567.  
  568. #ifdef PROTO
  569. static int getStr(FILE *inFile, char *token)
  570. #else
  571. static int getStr(inFile, token) 
  572. FILE *inFile;
  573. char *token;
  574. #endif
  575. {
  576.   int i;
  577.   int c;
  578.  
  579.   do 
  580.     {
  581.       while((c = fgetc(inFile)) != EOF && (c == ' ' || c == '\t'))
  582.          ;
  583.       if (c == '#')
  584.         while((c = fgetc(inFile)) != EOF && c != '\n')
  585.            ;
  586.     } while(c == '\n');
  587.  
  588.   for (i = 0; c != EOF && c != ' ' && c != '\t' && c != '\n'; i++)
  589.     {
  590.       token[i] = c;
  591.       c = fgetc(inFile);
  592.     }
  593.  
  594.   if (c != EOF)
  595.     {
  596.       ungetc(c, inFile);
  597.     }
  598.       
  599.   token[i] = '\0';
  600.   return i;
  601. } /* getStr */ 
  602.  
  603.  
  604. #ifdef PROTO
  605. static real *readValue(FILE *inFile, STF_spec *datasetSpec, ObjPtr dataset)
  606. #else
  607. static real *readValue(inFile, datasetSpec, dataset)
  608. FILE *inFile;
  609. STF_spec *datasetSpec;
  610. ObjPtr dataset;
  611. #endif
  612. /* return pointer to bounds of the dataset just read */
  613. /* return null if an error occurred */
  614. {
  615.  
  616.   real *bounds;
  617.   long *index;
  618.   int counter;
  619.   int nComponents;
  620.  
  621.   
  622.   /* determine number of components first */
  623.   if ( datasetSpec->fieldType == STF_SCALAR )
  624.     {
  625.       nComponents = 1;
  626.     }
  627.   else
  628.     {
  629.       nComponents = datasetSpec->nComponents;
  630.     }
  631.  
  632.  
  633.   /* alloc bounds array and init it to contain 0 values */
  634.   if ( !(bounds = (real *) Alloc( 2*nComponents*sizeof(real) )) )
  635.     {
  636.       FileFormatError("ReadSTFFile", "Unable to allocate memory");
  637.       return NULL;
  638.     }
  639.   for ( counter = 0; counter < 2*nComponents; counter++)
  640.      {
  641.        bounds[counter] = 0.0;
  642.      }
  643.  
  644.  
  645.   /* alloc index array and init it to contain 0 values */
  646.   if ( !(index = (long *) Alloc( datasetSpec->rank*sizeof(long) )) )
  647.     {
  648.       FileFormatError("ReadSTFFile", "Unable to allocate memory");
  649.       return 0;
  650.     }
  651.   for ( counter = 0; counter < datasetSpec->rank; counter++)
  652.      {
  653.        index[counter] = 0L;
  654.      }
  655.  
  656.  
  657.     /* prepare to read in the field values from file */
  658.     if ( !SetCurField(FIELD1, dataset) ) 
  659.       {
  660.         Free(index);
  661.         Free(bounds);
  662.         return NULL;
  663.       }
  664.  
  665.     /* adjust the index and read in the corresponding field value */ 
  666.     {
  667.       real curData;
  668.       char dataStr[32];
  669.       int n; /* counter */
  670.       
  671.       counter = 0;
  672.  
  673.       if ( datasetSpec->fieldType == STF_SCALAR ||
  674.                    datasetSpec->permutation == STF_NONINTERLACED )
  675.         {
  676.          for ( n = 0; n < nComponents; n++ )
  677.            {
  678.              counter = 0;
  679.              while(getStr(inFile, dataStr)) 
  680.              {
  681.                 int first = 1;
  682.                 int j;
  683.  
  684.                 if ( !ParseReal( &curData, dataStr ) )
  685.                   {
  686.                     FileFormatError("ReadSTFFile", "Invalid numeric Value");
  687.                     /* clean up and return an NULL ptr */
  688.                     Free(index);
  689.                     Free(bounds);
  690.                     return NULL;
  691.                   }
  692.  
  693.                 PutFieldComponent(FIELD1, n, index, (real) curData);
  694.  
  695.                 for (j = 0; first && j < datasetSpec->rank; j++)
  696.                    first = index[j] == 0;
  697.                 if (first)
  698.                   {
  699.                     bounds[2 * n] = curData;
  700.                     bounds[2 * n + 1] = curData;
  701.                   }
  702.                 else
  703.                   {
  704.                     /* looking for bounds */
  705.                     if ( curData < bounds[2*n] )
  706.                       {
  707.                         bounds[2*n] = curData;
  708.                       }
  709.                     else if ( curData > bounds[2*n+1] )
  710.                       {
  711.                         bounds[2*n+1] = curData;
  712.                       }
  713.                   }
  714.  
  715.                 /* adjust the index according to column or row major */
  716.                 if( datasetSpec->order == STF_COLUMNMAJOR)
  717.                   {
  718.                     for(counter = 0; counter < datasetSpec->rank; counter++)
  719.                       {
  720.                         index[counter]++;
  721.                         if (index[counter] == datasetSpec->dimensions[counter])
  722.                           {
  723.                            index[counter] = 0;
  724.                            continue;
  725.                           }
  726.                         break;
  727.                       }
  728.                     if (counter == datasetSpec->rank)
  729.                       {
  730.                         /* This signifies the nth component data is all read */
  731.                         break;
  732.                       }
  733.                   }
  734.                 else
  735.                   {
  736.                     for(counter = datasetSpec->rank - 1; counter > -1 ; counter--)
  737.                       {
  738.                         index[counter]++;
  739.                         if (index[counter] == datasetSpec->dimensions[counter])
  740.                           {
  741.                            index[counter] = 0;
  742.                            continue;
  743.                           }
  744.                         break;
  745.                       }
  746.                     if (counter == -1)
  747.                       {
  748.                         /* This signifies the nth component value is all read */
  749.                         break;
  750.                       }
  751.                   }
  752.  
  753.              } /* while */
  754.  
  755.            }  /* for ( n = 0; n < nComponents; n++ ) */
  756.  
  757.          /* this indicates that the data number is less than dimensions */
  758.          if ( ( datasetSpec->fieldType == STF_VECTOR && n != nComponents) ||
  759.               ( datasetSpec->order == STF_COLUMNMAJOR && counter != datasetSpec->rank )||
  760.               ( datasetSpec->order == STF_ROWMAJOR && counter != -1 ) )
  761.            {
  762.             FileFormatError("ReadSTFFile", "Inconsistent Data File, number of data values is less than that specified by dimensions");
  763.             Free(index);
  764.             Free(bounds);
  765.             return NULL;
  766.            }
  767.  
  768.         } /* if fieldType == STF_SCALAR || permutation == STF_NONINTERLACED */
  769.       else /* fieldType is VECTOR and permutation is INTERLACED */ 
  770.         {
  771.            n = 0;
  772.            counter = 0;
  773.            while(getStr(inFile, dataStr))
  774.            {
  775.              int first = 1;
  776.              int j;
  777.  
  778.               if ( !ParseReal( &curData, dataStr ) )
  779.                 {
  780.                   FileFormatError("ReadSTFFile", "Invalid numeric Value");
  781.                   /* clean up and return NULL */
  782.                   Free(index);
  783.                   Free(bounds);
  784.                   return NULL;
  785.                 }
  786.  
  787.               PutFieldComponent(FIELD1, n, index, (real) curData);
  788.  
  789.               for (j = 0; first && j < datasetSpec->rank; j++)
  790.                  first = index[j] == 0;
  791.               if (first)
  792.                 {
  793.                   bounds[2 * n] = curData;
  794.                   bounds[2 * n + 1] = curData;
  795.                 }
  796.               else
  797.                 {
  798.                   /* looking for bounds */
  799.                   if ( curData < bounds[2*n] )
  800.                     {
  801.                       bounds[2*n] = curData;
  802.                     }
  803.                   else if ( curData > bounds[2*n+1] )
  804.                     {
  805.                       bounds[2*n+1] = curData;
  806.                     }
  807.                  }
  808.  
  809.               n++;
  810.  
  811.               if ( n == nComponents ) /* time to update index */
  812.                 {
  813.                  /* adjust the index according to column or row major */
  814.                  if( datasetSpec->order == STF_COLUMNMAJOR)
  815.                    {
  816.                      for(counter = 0; counter < datasetSpec->rank; counter++)
  817.                        {
  818.                          index[counter]++;
  819.                          if (index[counter] == datasetSpec->dimensions[counter])
  820.                            {
  821.                             index[counter] = 0;
  822.                             continue;
  823.                            }
  824.                          break;
  825.                        }
  826.                      if (counter == datasetSpec->rank)
  827.                        {
  828.                          /* This signifies the field data is all read */
  829.                          break;
  830.                        }
  831.                    }
  832.                  else
  833.                    {
  834.                      for(counter = datasetSpec->rank - 1; counter > -1 ; counter--)
  835.                        {
  836.                          index[counter]++;
  837.                          if (index[counter] == datasetSpec->dimensions[counter])
  838.                            {
  839.                             index[counter] = 0;
  840.                             continue;
  841.                            }
  842.                          break;
  843.                        }
  844.                      if (counter == -1)
  845.                        {
  846.                          /* This signifies the field data is all read */
  847.                          break;
  848.                        }
  849.                    }
  850.  
  851.                  n = 0;
  852.  
  853.                 } /* if n == nComponents */
  854.  
  855.              } /* for fscanf */
  856.  
  857.          /* this indicates that the data number is less than dimensions */
  858.          if ( ( datasetSpec->order == STF_COLUMNMAJOR && counter != datasetSpec->rank )||
  859.               ( datasetSpec->order == STF_ROWMAJOR && counter != -1 ) )
  860.            {
  861.             FileFormatError("ReadSTFFile", "Inconsistent Data File, number of data values is less than that specified by dimensions");
  862.             Free(index);
  863.             Free(bounds);
  864.             return NULL;
  865.            }
  866.  
  867.         } /* else where fieldType is VECTOR and INTERLACED */
  868.  
  869.  
  870.     }
  871.  
  872.     Free(index);
  873.     return bounds;
  874.  
  875. }/* readValue */
  876.  
  877.  
  878. #ifdef PROTO
  879. static ObjPtr ReadSTFFile(char *name)
  880. #else
  881. static ObjPtr ReadSTFFile(name)
  882. char *name;
  883. #endif
  884. {
  885.     int datasetType; /* DATA, GRID, END or ERROR */
  886.     STF_spec datasetSpec;
  887.     int specFlag = 0;
  888.  
  889.     FILE *inFile;
  890.  
  891.     real *gridBounds = NULL;
  892.     real *dataMinmaxes = NULL;
  893.  
  894.     int isGrid = 0, isData = 0;
  895.     int gridComponents; /* the number of components for the grid if any */
  896.     int counter;
  897.     int more = 1;
  898.  
  899.     ObjPtr dataForm;
  900.     ObjPtr gridDataset;
  901.     ObjPtr retVal;
  902.  
  903.     inFile = fopen(name, "r");
  904.     if ( !inFile )
  905.       {
  906.         Error("ReadSTFFile", OPENFILEERROR, name);
  907.         return NULLOBJ;
  908.       }
  909.      
  910.  
  911.     /* set default dataset specs */
  912.     datasetSpec.rank = 0;
  913.     datasetSpec.nComponents = 0;
  914.     datasetSpec.order = STF_COLUMNMAJOR;
  915.     datasetSpec.permutation = STF_INTERLACED;
  916.     strcpy(datasetSpec.datasetName, name);
  917.  
  918.  
  919.  
  920.     while( more ) 
  921.       {
  922.         /* parse header lines */
  923.         datasetType = parseHeader(inFile, &specFlag, &datasetSpec);
  924.  
  925.         switch ( datasetType )
  926.           {
  927.             case STF_ERROR: 
  928.                FileFormatError("Error in dataset", datasetSpec.datasetName);
  929.                more = 0;
  930.                break;
  931.             case STF_GRID:
  932.                /* if isGrid then error, double grid */
  933.                if ( isGrid )
  934.                  {
  935.                    FileFormatError("ReadSTFFile", "Multiple Grid Spec");
  936.                    FileFormatError("Error in dataset", datasetSpec.datasetName);
  937.                    more = 0;
  938.                    break;
  939.                  }
  940.  
  941.                gridComponents = datasetSpec.nComponents;
  942.  
  943.                /* create the grid dataset */
  944.                gridDataset = NewStructuredDataset(datasetSpec.datasetName, datasetSpec.rank, datasetSpec.dimensions, gridComponents);
  945.  
  946.                /* something wrong with grid value format */
  947.                if ( !(gridBounds = readValue( inFile, &datasetSpec, gridDataset)))
  948.                  {
  949.                    FileFormatError("Error in dataset", datasetSpec.datasetName);
  950.                    more = 0;
  951.                    break;
  952.                  }
  953.  
  954.                /* reset specFlag */
  955.                specFlag &= (~STF_ISTYPE & ~STF_ISPERMUTATION & ~STF_ISORDER);
  956.                isGrid = 1;
  957.                break;
  958.             case STF_DATA:
  959.                /* if isData then error, double data */
  960.                if ( isData )
  961.                  {
  962.                    FileFormatError("ReadSTFFile", "Multiple Data Spec");
  963.                    FileFormatError("Error in dataset", datasetSpec.datasetName);
  964.                    more = 0;
  965.                    break;
  966.                  }
  967.  
  968.                /* create the dataset */
  969.                retVal = NewStructuredDataset(datasetSpec.datasetName, datasetSpec.rank, datasetSpec.dimensions, datasetSpec.nComponents);
  970.  
  971.                /* something wrong with data value format */
  972.                if ( !(dataMinmaxes = readValue( inFile, &datasetSpec, retVal)) )
  973.                  {
  974.                    FileFormatError("Error in dataset", datasetSpec.datasetName);
  975.                    more = 0;
  976.                    break;
  977.                  }
  978.                else
  979.                  {
  980.                    /* has no use with dataMinmaxes */
  981.                    Free(dataMinmaxes);
  982.                  }
  983.  
  984.                /* reset specFlag */
  985.                specFlag &= (~STF_ISTYPE & ~STF_ISPERMUTATION & ~STF_ISORDER);
  986.                isData = 1;
  987.                break;
  988.             case STF_END:
  989.                /* if isData && isGrid then use curvilinear dataForm */
  990.                /* if isData && !isGrid then regular dataForm */
  991.                /* for this case, the bounds are required, error if no */
  992.                /* if !isData && isGrid the error */
  993.                /* if !isData && !isGrid then means end of file */
  994.                        
  995.                if ( !isData && !isGrid )
  996.                  /* end of dataset file reached */
  997.                  {
  998.                    more = 0;
  999.                    break;
  1000.                  }
  1001.  
  1002.                if ( !isData && isGrid )
  1003.                  {
  1004.                    FileFormatError("ReadSTFFile", "NO DATA");
  1005.                    FileFormatError("Error in dataset", datasetSpec.datasetName);
  1006.                    more = 0;
  1007.                    break;
  1008.                  }
  1009.  
  1010.                if ( specFlag & STF_ISBOUNDS )
  1011.                  {
  1012.                    int nBounds; /* number of bounds */
  1013.  
  1014.                    /* get rid of calculated grid bounds if any and 
  1015.                       use the bounds provided by the user */
  1016.                    /* determine the number of bounds */
  1017.                    if ( isGrid )
  1018.                      {
  1019.                        Free(gridBounds);
  1020.                        gridBounds = NULL;
  1021.                        nBounds = 2 * gridComponents;
  1022.                      }
  1023.                    else
  1024.                      {
  1025.                        nBounds = 2 * datasetSpec.rank;
  1026.                      }
  1027.  
  1028.                    /* test if the user provide correct bounds */
  1029.                    if ( nBounds != datasetSpec.nBounds )
  1030.                      {
  1031.                        FileFormatError("ReadSTFFile", "BOUNDS format error");
  1032.                        FileFormatError("Error in dataset", datasetSpec.datasetName);
  1033.                        more = 0;
  1034.                        break;
  1035.                      }
  1036.  
  1037.                    gridBounds = datasetSpec.bounds;
  1038.                  }
  1039.  
  1040.  
  1041.                if ( isData && isGrid )
  1042.                  /* for this case, bounds are optional */
  1043.                  {
  1044.                    dataForm = NewCurvilinearDataForm(datasetSpec.datasetName, datasetSpec.rank, datasetSpec.dimensions, gridBounds, gridDataset);
  1045.                  }
  1046.                else
  1047.                  {
  1048.                    /* isData && !isGrid */      
  1049.                    /* for this case bounds are required */
  1050.                    if ( !(specFlag & STF_ISBOUNDS) )
  1051.                      {
  1052.                        FileFormatError("ReadSTFFile", "BOUNDS required for dataset with no grid");
  1053.                        FileFormatError("Error in dataset", datasetSpec.datasetName);
  1054.                        more = 0;
  1055.                        break;
  1056.                      }
  1057.                    dataForm = NewRegularDataForm(datasetSpec.datasetName, datasetSpec.rank, datasetSpec.dimensions, gridBounds);
  1058.                  }
  1059.  
  1060.                SetDatasetForm(retVal, dataForm);
  1061.                if ( specFlag & STF_ISTIME )
  1062.                  {
  1063.                    RegisterTimeDataset(retVal, datasetSpec.time);
  1064.                  }
  1065.                else
  1066.                  {
  1067.                    RegisterDataset(retVal);
  1068.                  }
  1069.  
  1070.                isGrid = 0;
  1071.                isData = 0;
  1072.                specFlag = 0;
  1073.                Free(datasetSpec.dimensions);
  1074.                Free(gridBounds);
  1075.                gridBounds = NULL;
  1076.  
  1077.                /* set default dataset specs */
  1078.                datasetSpec.rank = 0;
  1079.                datasetSpec.nComponents = 0;
  1080.                datasetSpec.order = STF_COLUMNMAJOR;
  1081.            datasetSpec.permutation = STF_INTERLACED;
  1082.  
  1083.                break;
  1084.             default:     
  1085.                /* dummy case */
  1086.                break;
  1087.           } /* switch */
  1088.       } /* while */
  1089.  
  1090.    if ( specFlag & STF_ISDIMS )
  1091.      {
  1092.        Free(datasetSpec.dimensions);
  1093.      }
  1094.  
  1095.    if ( specFlag & STF_ISBOUNDS )
  1096.      {
  1097.        Free(datasetSpec.bounds);
  1098.      }
  1099.    if ( isGrid )
  1100.      {
  1101.        Free(gridBounds);
  1102.      }
  1103.    fclose(inFile);
  1104.    return NULLOBJ;
  1105. }
  1106.  
  1107. /*******************************************************************************/
  1108. /********************************P3D READER ************************************/
  1109. /*******************************************************************************/
  1110. /* Tzong-Yow Hwu */
  1111. /* 10/26/1992 */
  1112. /* plot3d file reader */
  1113. /* 9/7/1993, added SGI binary format */
  1114.  
  1115. /* file types */
  1116. #define P3D_GRID        0
  1117. #define P3D_FUNC        1
  1118. #define P3D_Q           2
  1119.  
  1120. /* file format types, default on IRIS 4D is unformatted and on UNIX is binary */
  1121. /* for both grid and Q files */
  1122. #define P3D_DAT         0  /* default */ 
  1123. #define P3D_BIN         1 
  1124. #define P3D_FMT         2
  1125.  
  1126. /* data permutation method, default is the whole */
  1127. /* PLANES has no effect on reading other than 3D dataset */
  1128. #define P3D_WHOLE       0   /* default */   /* the whole dataset at once */
  1129. #define P3D_PLANE       1   /* one plane at a time */
  1130.  
  1131. /* Control over default checking of Q data, default is to check */
  1132. /* for Q files only */
  1133. #define P3D_CHECK       0  /* default */
  1134. #define P3D_NOCHECK     1
  1135.  
  1136.  
  1137. /* Jacobian control */
  1138. /* for Q files only */
  1139. #define P3D_NOJACOBIAN  0 /* default */
  1140. #define P3D_JACOBIAN    1 
  1141.  
  1142. /* IBLANK control */
  1143. /* for grid files only */
  1144. #define P3D_NOBLANK     0
  1145. #define P3D_BLANK       1
  1146.  
  1147. #define MAXNAMELEN 1800 
  1148.  
  1149. typedef struct P3D_grid
  1150. {
  1151.   char name[MAXNAMELEN];
  1152.   int fileType;
  1153.   int permutation;
  1154.   int blank;
  1155. } p3d_grid;
  1156.  
  1157. typedef struct P3D_q
  1158. {
  1159.   char name[MAXNAMELEN];
  1160.   int isMDataset;
  1161.   int fileType;
  1162.   int permutation;
  1163.   int check;
  1164.   int jacobian;
  1165. } p3d_q;
  1166.  
  1167. typedef struct P3D_func
  1168. {
  1169.   char name[MAXNAMELEN];
  1170.   int fileType;
  1171.   int permutation;
  1172. } p3d_func;
  1173.  
  1174. /* a p3d set consists of a grid spec and several Q and/or func spec's */
  1175. typedef struct P3D_set
  1176.    {
  1177.      int ndim;
  1178.      int isMGrid;    /* 0 for no and 1 for yes */
  1179.      int isGrid;     /* is the grid present */
  1180.      p3d_grid *grid; /* pointer to the grid spec */
  1181.      int nQ;         /* number of Q files for this p3d set */
  1182.      p3d_q *qlist;   /* pointer to an array of Q spec */
  1183.      int nFunc;      /* number of func files for this p3d set */
  1184.      p3d_func *funclist; /* pointer to an array of func spec */
  1185.    } p3d_set;
  1186.  
  1187. /*
  1188.    1. getToken read to next occurrence of slash '/' or eof, and pass
  1189.       the string back, return 0 if end of file is encountered, -1 
  1190.       if error format, otherwise the number of chars in token.
  1191.    2. only process READ commands 
  1192.    3. newline without preceding '-' is treated as command terminator
  1193. */
  1194. #ifdef PROTO
  1195. static int getToken(FILE *inFile, char *token)
  1196. #else
  1197. static int getToken(inFile, token) 
  1198. FILE *inFile;
  1199. char *token;
  1200. #endif
  1201. {
  1202.   int i;
  1203.   int c;
  1204.   int inQuotes = 0;
  1205.  
  1206.   for ( i = 0; ( c = fgetc(inFile) ) != EOF ; i++ ) 
  1207.      {
  1208.        if (c == ' ' || c == '\t')
  1209.          {
  1210.           i--;
  1211.           continue;
  1212.          }
  1213.        else if (c == '\n')
  1214.          {
  1215.            if (i > 0 && token[i - 1] == '-')
  1216.              /* if continuation on next line */
  1217.              {
  1218.                i -= 2;
  1219.                continue;
  1220.              }
  1221.            else if (inQuotes)
  1222.              {
  1223.                fprintf(stderr, "File format error, unbalanced quote\n");
  1224.                return -1;
  1225.              }
  1226.            else if (i == 0)
  1227.              /* skip blank lines */
  1228.              {
  1229.                i--;
  1230.                continue;
  1231.              }
  1232.            else
  1233.              {
  1234.                break;
  1235.              }
  1236.          }
  1237.        else if (c == '\"')
  1238.          {
  1239.            /* flip inQuotes flag */
  1240.            inQuotes = !inQuotes;
  1241.            /* discard the '"' char */
  1242.            i--;
  1243.            continue;
  1244.          }
  1245.        else if (c == '/')
  1246.          {
  1247.            if (!inQuotes)
  1248.              {
  1249.                break;
  1250.              }
  1251.            /* otherwise, fall through to get the char */
  1252.          }
  1253.        token[i] = c;
  1254.      }
  1255.   
  1256.   if ( i == 0 )
  1257.     if ( c == EOF )
  1258.       {
  1259.         /* end of file reached */
  1260.         return 0;
  1261.       }
  1262.     else
  1263.       {
  1264.         /* error file formate */
  1265.         /* call ERROR routine */
  1266.         fprintf(stderr, "File format error\n");
  1267.         return -1;
  1268.       }
  1269.  
  1270.   token[i] = '\0';
  1271.   return i;
  1272. } /* getToken */ 
  1273.  
  1274. /* getvalue return the number of values repeated in the form
  1275.    of num*val and returns 0 if error or EOF encountered */
  1276. #ifdef PROTO
  1277. static int getValue(FILE *inFile, double *val)
  1278. #else
  1279. static int getValue(inFile, val)
  1280. FILE *inFile;
  1281. double *val;
  1282. #endif
  1283. {
  1284.   int c;
  1285.   int count;
  1286.   char s[180];
  1287.   char *end, *newstr;
  1288.  
  1289.   /* skip blanks and commas */
  1290.   while(isspace(c = fgetc(inFile)) || c == ',')
  1291.      ;
  1292.   if (c != EOF)
  1293.     ungetc(c, inFile);
  1294.   else
  1295.     return EOF;
  1296.  
  1297.   /* get the numerical string */
  1298.   if (fscanf(inFile, "%s", s) != 1)
  1299.     {
  1300.       return 0;
  1301.     }
  1302.   *val = strtod(s, &end);
  1303.   if (end == s)
  1304.     {
  1305.       /* nothing is converted */
  1306.       return 0;
  1307.     }
  1308.   if (end[0] == '\0' || (end[0] == ',' && end[1] == '\0'))
  1309.     {
  1310.       return 1;
  1311.     }
  1312.   else if (*end == '*')
  1313.     {
  1314.       count = (int) *val;
  1315.       newstr = end + 1;
  1316.       *val = strtod(newstr, &end);
  1317.       if (end == newstr || !(*end == '\0' || (*end = ',' && *(end + 1) == '\0')))
  1318.         {
  1319.           /* nothing is converted || format error */
  1320.           return 0;
  1321.         }
  1322.       return count;
  1323.     }
  1324.   else
  1325.     {
  1326.       return 0;
  1327.     }
  1328. } /* getValue */
  1329.  
  1330. /* strncmp2 compares two string for n chars, case insensitive */
  1331. #ifdef PROTO
  1332. static int strncmp2(char *s, char *t, int n)
  1333. #else
  1334. static int strncmp2(s, t, n)
  1335. char *s;
  1336. char *t;
  1337. int n;
  1338. #endif
  1339. {
  1340.   int i;
  1341.  
  1342.   for (i = 0; i < n && toupper(s[i]) == toupper(t[i]); i++)
  1343.      {
  1344.        if (s[i] == '\0')
  1345.          return 0;
  1346.      }
  1347.  
  1348.   if (i == n)
  1349.     return 0;
  1350.   else
  1351.     return (s[i] - t[i]);
  1352. }
  1353.  
  1354. #ifdef PROTO
  1355. static void cleanup(int nSets, p3d_set *sets)
  1356. #else
  1357. static void cleanup(nSets, sets) 
  1358. int nSets;
  1359. p3d_set *sets;
  1360. #endif
  1361. /* cleanup free the memory alloced in the program */
  1362. {
  1363.   int i;
  1364.  
  1365.   if (!sets)
  1366.     return;
  1367.  
  1368.   for (i = 0; i < nSets; i++)
  1369.      {
  1370.        if ((sets + i)->nQ)
  1371.           {
  1372.             Free((sets + i)->qlist);
  1373.           }
  1374.        if ((sets + i)->nFunc)
  1375.           {
  1376.             Free((sets + i)->funclist);
  1377.           }
  1378.        if ((sets + i)->isGrid)
  1379.          {
  1380.             Free((sets + i)->grid);
  1381.          }
  1382.      }
  1383.  
  1384.   Free(sets);
  1385.  
  1386.   return;
  1387. } /* cleanup */ 
  1388.  
  1389. #ifdef PROTO
  1390. static void compress(int *nSets, p3d_set **sets)
  1391. #else
  1392. static void compress(nSets, sets) 
  1393. int *nSets;
  1394. p3d_set **sets;
  1395. #endif
  1396. {
  1397.   int i, j;
  1398.  
  1399.   /* now, sort the sets according to their grid filename first */
  1400.   for (i = 1; i < *nSets; i++)
  1401.     {
  1402.      p3d_set temp;
  1403.  
  1404.      temp = (*sets)[i];
  1405.      /* filename is case sensitive, so use strcmp not strcmp2 */
  1406.      for (j = i - 1; j >= 0 && strcmp((*sets)[j].grid->name,
  1407.                                        temp.grid->name) > 0; j--)
  1408.         {
  1409.           (*sets)[j + 1] = (*sets)[j];
  1410.         }
  1411.      (*sets)[j + 1] = temp;
  1412.     }
  1413.  
  1414.   /* second compress the sets, so that if any of these sets use the
  1415.      same grid file, making them as the same set.
  1416.   */
  1417.   {
  1418.     p3d_set *last = *sets;
  1419.     p3d_set *temp;
  1420.  
  1421.     for (i = 1, j = 1; i < *nSets; i++)
  1422.       {
  1423.         p3d_set *curset = *sets + i;
  1424.  
  1425.         if (!strcmp(curset->grid->name, last->grid->name))
  1426.           {
  1427.             /* combine this set into last set */
  1428.             if (curset->nQ > 0)
  1429.               {
  1430.                 int k;
  1431.  
  1432.                 if (!(last->qlist = (p3d_q *) Realloc(last->qlist, 
  1433.                                   (last->nQ + curset->nQ) * sizeof(p3d_q))))
  1434.                   {
  1435.                     fprintf(stderr, "Unable to alloc memory, aborted\n");
  1436.                     return;
  1437.                   }
  1438.                 for (k = 0; k < curset->nQ; k++)
  1439.                   {
  1440.                     last->qlist[last->nQ + k] = curset->qlist[k];
  1441.                   }
  1442.                 Free(curset->qlist);
  1443.                 curset->qlist = NULL;
  1444.                 last->nQ += curset->nQ;
  1445.               }
  1446.  
  1447.             if (curset->nFunc > 0)
  1448.               {
  1449.                 int k;
  1450.  
  1451.                 if (!(last->funclist = (p3d_func *) Realloc(last->funclist, 
  1452.                             (last->nFunc + curset->nFunc) * sizeof(p3d_func))))
  1453.                   {
  1454.                     fprintf(stderr, "Unable to alloc memory, aborted\n");
  1455.                     return;
  1456.                   }
  1457.                 for (k = 0; k < curset->nFunc; k++)
  1458.                   {
  1459.                     last->funclist[last->nFunc + k] = curset->funclist[k];
  1460.                   }
  1461.                 Free(curset->funclist);
  1462.                 curset->funclist = NULL;
  1463.                 last->nFunc += curset->nFunc;
  1464.               }
  1465.             
  1466.             /* mark this set has been combined */
  1467.             Free(curset->grid);
  1468.             curset->grid = NULL;
  1469.           }
  1470.         else
  1471.           {
  1472.             last[j++] = *curset; 
  1473.           }
  1474.       } /* for */
  1475.       
  1476.  
  1477.     /* reset sets and nSets */
  1478.     /* cannot use realloc for *sets, since if it fails, then the q and func lists 
  1479.        is unaccessible and becomes garbages occupying the memory */
  1480.     if (!(temp = (p3d_set *) Alloc(j * sizeof(p3d_set))))
  1481.       {
  1482.         fprintf(stderr, "Unable to alloc memory, aborted\n");
  1483.         return;
  1484.       }
  1485.     for (i = 0; i < j; i++)
  1486.       {
  1487.         temp[i] = (*sets)[i];
  1488.       }
  1489.     Free(*sets); 
  1490.     *sets = temp;
  1491.     *nSets = j;
  1492.   }
  1493.  
  1494. #ifdef DEBUG
  1495.  
  1496.   for (i = 0; i < *nSets; i++)
  1497.     {
  1498.  
  1499.       p3d_set *curset = *sets + i;
  1500.       int isGrid = curset->isGrid;
  1501.       p3d_grid *grid = curset->grid;
  1502.       int nQ = curset->nQ;
  1503.       p3d_q *qlist = curset->qlist;
  1504.       int nFunc = curset->nFunc;
  1505.       p3d_func *funclist = curset->funclist;
  1506.       int j;
  1507.  
  1508.       printf("READ");
  1509.       /* ndim */
  1510.       if (curset->ndim == 1)
  1511.         {
  1512.           printf("/1D"); 
  1513.         }
  1514.       else if (curset->ndim == 2)
  1515.         {
  1516.           printf("/2D");
  1517.         }
  1518.       else if (curset->ndim == 3)
  1519.         {
  1520.           printf("/3D");
  1521.         }
  1522.     
  1523.       /* MGRID */
  1524.       if (curset->isMGrid)
  1525.         {
  1526.           printf("/MGRID");
  1527.         }
  1528.  
  1529.       /* process the XYZ file */
  1530.       if (isGrid)
  1531.         {
  1532.           /* file format */
  1533.           if (grid->fileType == P3D_DAT)
  1534.             {
  1535.               printf("/UNFORMATTED");
  1536.             }
  1537.           else if (grid->fileType == P3D_BIN)
  1538.             {
  1539.               printf("/BINARY");
  1540.             }
  1541.           else if (grid->fileType == P3D_FMT)
  1542.             {
  1543.               printf("/FORMATTED");
  1544.             }
  1545.  
  1546.           /* permutation */
  1547.           if (grid->permutation == P3D_WHOLE)
  1548.             {
  1549.               printf("/WHOLE");
  1550.             }
  1551.           else if (grid->permutation == P3D_PLANE)
  1552.             {
  1553.               printf("/PLANES");
  1554.             }
  1555.  
  1556.           /* BLANK */
  1557.           if (grid->permutation == P3D_BLANK)
  1558.             {
  1559.               printf("/BLANK");
  1560.             }
  1561.           else if (grid->permutation == P3D_NOBLANK)
  1562.             {
  1563.               printf("/NOBLANK");
  1564.             }
  1565.  
  1566.           printf("/XYZ=%s", grid->name);
  1567.         }
  1568.       else
  1569.          /* should not be the case */
  1570.         {
  1571.           fprintf(stderr, "XYZ file missing\n");
  1572.           continue;  /* do nothing for this set */
  1573.         }
  1574.  
  1575.       /* process the Q files */
  1576.       for (j = 0; j < nQ; j++)
  1577.         {
  1578.           /* MDATASET */
  1579.           if (qlist[j].isMDataset)
  1580.             {
  1581.               printf("/MDATASET");
  1582.             }
  1583.  
  1584.           /* file format */
  1585.           if (qlist[j].fileType == P3D_DAT)
  1586.             {
  1587.               printf("/UNFORMATTED");
  1588.             }
  1589.           else if (qlist[j].fileType == P3D_BIN)
  1590.             {
  1591.               printf("/BINARY");
  1592.             }
  1593.           else if (qlist[j].fileType == P3D_FMT)
  1594.             {
  1595.               printf("/FORMATTED");
  1596.             }
  1597.  
  1598.           /* permutation */
  1599.           if (qlist[j].permutation == P3D_WHOLE)
  1600.             {
  1601.               printf("/WHOLE");
  1602.             }
  1603.           else if (qlist[j].permutation == P3D_PLANE)
  1604.             {
  1605.               printf("/PLANES");
  1606.             }
  1607.  
  1608.           /* CHECK */
  1609.           if (qlist[j].check)
  1610.             {
  1611.               printf("/CHECK");
  1612.             }
  1613.           else
  1614.             {
  1615.               printf("/NOCHECK");
  1616.             }
  1617.  
  1618.           /* JACOBIAN */
  1619.           if (qlist[j].jacobian)
  1620.             {
  1621.               printf("/JACOBIAN");
  1622.               printf("\nsorry, jacobian is not implemented\n");
  1623.             }
  1624.           else
  1625.             {
  1626.               printf("/NOJACOBIAN");
  1627.             }
  1628.  
  1629.           printf("/Q=%s", qlist[j].name);
  1630.         }
  1631.  
  1632.       /* process the FUNCTION files */
  1633.       for (j = 0; j < nFunc; j++)
  1634.         {
  1635.           /* file format */
  1636.           if (funclist[j].fileType == P3D_DAT)
  1637.             {
  1638.               printf("/UNFORMATTED");
  1639.             }
  1640.           else if (funclist[j].fileType == P3D_BIN)
  1641.             {
  1642.               printf("/BINARY");
  1643.             }
  1644.           else if (funclist[j].fileType == P3D_FMT)
  1645.             {
  1646.               printf("/FORMATTED");
  1647.             }
  1648.  
  1649.           /* permutation */
  1650.           if (funclist[j].permutation == P3D_WHOLE)
  1651.             {
  1652.               printf("/WHOLE");
  1653.             }
  1654.           else if (funclist[j].permutation == P3D_PLANE)
  1655.             {
  1656.               printf("/PLANES");
  1657.             }
  1658.  
  1659.           printf("/FUNCTION=%s", funclist[j].name);
  1660.         }
  1661.  
  1662.       printf("\n");
  1663.     } /* outter for */
  1664.  
  1665. #endif   /* for ifdef DEBUG */
  1666.  
  1667.   return;
  1668. } /* compress */ 
  1669.  
  1670.  
  1671. /* return 0 on failure */
  1672. #ifdef PROTO
  1673. static int parseP3DMeta(char *p3dname, int *nSets, p3d_set **sets)
  1674. #else
  1675. static int parseP3DMeta(p3dname, nSets, sets)
  1676. char *p3dname;
  1677. int *nSets;
  1678. p3d_set **sets;
  1679. #endif
  1680. {
  1681.  
  1682.   FILE *inFile;
  1683.   char token[MAXNAMELEN];
  1684.   int isNdim = 0, ndim = 3;  /* default is 3D */
  1685.   int isMDataset = 0;  /* time dependent datasets */
  1686.   int isMGrid = 0;    /* multiple grid files */  /* How do deal with it in Scian? */
  1687.   int isFileType = 0, fileType = P3D_DAT; /* file type, default is UNFORMATTED */
  1688.   int isPermutation = 0, permutation = P3D_WHOLE; /* data permutation way */
  1689.   int isCheck = 0, check = P3D_NOCHECK;  /* check control for Q file */
  1690.   int isJacobian = 0, jacobian = P3D_NOJACOBIAN;  /* jacobian control for grid file*/
  1691.   int isBlank = 0, blank = P3D_NOBLANK;  /* blank control for grid file*/
  1692.   int status;
  1693.   int valid = 1; /* if this file format is valid */
  1694.  
  1695.  
  1696.   if ( !(inFile = fopen(p3dname, "r")) )
  1697.     {
  1698.       fprintf(stderr, "ReadP3DFile: %s\n", p3dname);
  1699.       return 0;
  1700.     }
  1701.  
  1702.   /* the first thing the .p3d file should be READ */
  1703.   if ((status = getToken(inFile, token)) <= 0 || strcmp2(token, "READ"))
  1704.     { 
  1705.       fprintf(stderr, "ReadP3DFile: %s, file format error\n", p3dname);
  1706.       fclose(inFile);
  1707.       return 0;
  1708.     }
  1709.  
  1710.   /* now since we know that this is plot 3d file, lets process it */
  1711.   /* alloc memory for the first p3d set */
  1712.   if (!(*sets = (p3d_set *) Alloc(sizeof(p3d_set))))
  1713.     {
  1714.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  1715.       fclose(inFile);
  1716.       return 0;
  1717.     }
  1718.   (*sets)->ndim = 3;
  1719.   (*sets)->isMGrid = 0;
  1720.   (*sets)->isGrid = 0;
  1721.   (*sets)->grid = NULL;
  1722.   (*sets)->nQ = 0;
  1723.   (*sets)->qlist = NULL;
  1724.   (*sets)->nFunc = 0;
  1725.   (*sets)->funclist = NULL;
  1726.   *nSets = 1;
  1727.  
  1728.   while ( (status = getToken(inFile, token)) > 0 )
  1729.     {
  1730.        /* for dimensions */
  1731.        if (!strcmp2(token, "1D") || !strcmp2(token, "1") ||
  1732.            !strcmp2(token, "2D") || !strcmp2(token, "2") ||
  1733.            !strcmp2(token, "3D") || !strcmp2(token, "3"))
  1734.          {
  1735.            p3d_set *curset = *sets + *nSets - 1;
  1736.  
  1737.            /* all number of dimensions for the same set must be consistent */
  1738.            if (isNdim && ndim != token[0] - '0')
  1739.              {
  1740.                fprintf(stderr, "ReadP3DFile: %s, inconsistent dim spec.\n", p3dname);
  1741.                valid = 0;
  1742.                break; /* get out of the while loop */ 
  1743.              }
  1744.  
  1745.            isNdim = 1;
  1746.            ndim = token[0] - '0';
  1747.            continue;
  1748.          }
  1749.  
  1750.        /* for MDATASET */
  1751.        if (!strcmp2(token, "MD") || !strcmp2(token, "MDA") ||
  1752.            !strcmp2(token, "MDAT") || !strcmp2(token, "MDATA") ||
  1753.            !strcmp2(token, "MDATAS") || !strcmp2(token, "MDATASE") ||
  1754.            !strcmp2(token, "MDATASET"))
  1755.          {
  1756.            isMDataset = 1;
  1757.            continue;
  1758.          }
  1759.  
  1760.        /* for MGRID */
  1761.        if (!strcmp2(token, "MG") || !strcmp2(token, "MGR") ||
  1762.            !strcmp2(token, "MGRI") || !strcmp2(token, "MGRID"))
  1763.          {
  1764.            /* all files for this set must be specified as mgrid */
  1765.            isMGrid = 1;
  1766.            continue;
  1767.          }
  1768.  
  1769.        /* for file format */
  1770.        if (!strcmp2(token, "FO") || !strcmp2(token, "FOR") ||
  1771.            !strcmp2(token, "FORM") || !strcmp2(token, "FORMA") || 
  1772.            !strcmp2(token, "FORMAT") || !strcmp2(token, "FORMATT") || 
  1773.            !strcmp2(token, "FORMATTE") || !strcmp2(token, "FORMATTED"))
  1774.  
  1775.          {
  1776.            if (isFileType && fileType != P3D_FMT)
  1777.              {
  1778.                 fprintf(stderr, "ReadP3DFile: %s, double file format spec.\n", p3dname);
  1779.                 valid = 0;
  1780.                 break; /* get out of the while loop */ 
  1781.              }
  1782.            isFileType = 1;
  1783.            fileType = P3D_FMT;
  1784.            continue;
  1785.          }
  1786.        if (!strcmp2(token, "U") || !strcmp2(token, "UN") ||
  1787.            !strcmp2(token, "UNF") || !strcmp2(token, "UNFO") ||
  1788.            !strcmp2(token, "UNFOR") || !strcmp2(token, "UNFORM") ||
  1789.            !strcmp2(token, "UNFORMA") || !strcmp2(token, "UNFORMAT") ||
  1790.            !strcmp2(token, "UNFORMATT") || !strcmp2(token, "UNFORMATTE") ||
  1791.            !strcmp2(token, "UNFORMATTED"))
  1792.          {
  1793.            if (isFileType && fileType != P3D_DAT)
  1794.              {
  1795.                 fprintf(stderr, "ReadP3DFile: %s, double file format spec.\n", p3dname);
  1796.                 valid = 0;
  1797.                 break; /* get out of the while loop */ 
  1798.              }
  1799.            isFileType = 1;
  1800.            fileType = P3D_DAT;
  1801.            continue;
  1802.          }
  1803.        if (!strcmp2(token, "BI") || !strcmp2(token, "BIN") ||
  1804.            !strcmp2(token, "BINA") || !strcmp2(token, "BINAR") ||
  1805.            !strcmp2(token, "BINARY"))
  1806.          {
  1807.            if (isFileType && fileType != P3D_BIN)
  1808.              {
  1809.                 fprintf(stderr, "ReadP3DFile: %s, double file format spec.\n", p3dname);
  1810.                 valid = 0;
  1811.                 break; /* get out of the while loop */ 
  1812.              }
  1813.            isFileType = 1;
  1814.            fileType = P3D_BIN;
  1815.            continue;
  1816.          }
  1817.  
  1818.  
  1819.        /* for data permutation */
  1820.        if (!strcmp2(token, "P") || !strcmp2(token, "PL") ||
  1821.            !strcmp2(token, "PLA") || !strcmp2(token, "PLAN") ||
  1822.            !strcmp2(token, "PLANE") || !strcmp2(token, "PLANES"))
  1823.          {
  1824.            if (isPermutation && permutation != P3D_PLANE)
  1825.              {
  1826.                 fprintf(stderr, "ReadP3DFile: %s, double permutation spec.\n", p3dname);
  1827.                 valid = 0;
  1828.                 break; /* get out of the while loop */ 
  1829.              }
  1830.            isPermutation = 1;
  1831.            permutation = P3D_PLANE;
  1832.            continue;
  1833.          }
  1834.        if (!strcmp2(token, "W") || !strcmp2(token, "WH") ||
  1835.            !strcmp2(token, "WHO") || !strcmp2(token, "WHOL") ||
  1836.            !strcmp2(token, "WHOLE"))
  1837.          {
  1838.            if (isPermutation && permutation != P3D_WHOLE)
  1839.              {
  1840.                 fprintf(stderr, "ReadP3DFile: %s, double permutation spec.\n", p3dname);
  1841.                 valid = 0;
  1842.                 break; /* get out of the while loop */ 
  1843.              }
  1844.            isPermutation = 1;
  1845.            permutation = P3D_WHOLE;
  1846.            continue;
  1847.          }
  1848.  
  1849.        /* for CHECK or not */
  1850.        if (!strcmp2(token, "C") || !strcmp2(token, "CH") ||
  1851.            !strcmp2(token, "CHE") || !strcmp2(token, "CHEC") ||
  1852.            !strcmp2(token, "CHECK"))
  1853.          {
  1854.            if (isCheck && check != P3D_CHECK)
  1855.              {
  1856.                 fprintf(stderr, "ReadP3DFile: %s, double check spec.\n", p3dname);
  1857.                 valid = 0;
  1858.                 break; /* get out of the while loop */ 
  1859.              }
  1860.            isCheck = 1;
  1861.            check = P3D_CHECK;
  1862.            continue;
  1863.          }
  1864.        if (!strcmp2(token, "NOC") || !strcmp2(token, "NOCH") ||
  1865.            !strcmp2(token, "NOCHE") || !strcmp2(token, "NOCHEC") ||
  1866.            !strcmp2(token, "NOCHECK"))
  1867.          {
  1868.            if (isCheck && check != P3D_NOCHECK)
  1869.              {
  1870.                 fprintf(stderr, "ReadP3DFile: %s, double check spec.\n", p3dname);
  1871.                 valid = 0;
  1872.                 break; /* get out of the while loop */ 
  1873.              }
  1874.            isCheck = 1;
  1875.            check = P3D_NOCHECK;
  1876.            continue;
  1877.          }
  1878.  
  1879.        /* for Jacobian spec */
  1880.            /* jacobian is not implemented */
  1881.        if (!strcmp2(token, "J") || !strcmp2(token, "JA") ||
  1882.            !strcmp2(token, "JAC") || !strcmp2(token, "JACO") ||
  1883.            !strcmp2(token, "JACOB") || !strcmp2(token, "JACOBI") ||
  1884.            !strcmp2(token, "JACOBIA") || !strcmp2(token, "JACOBIAN"))
  1885.          {
  1886.            /* 
  1887.            isJacobian = 1;
  1888.            check = P3D_JACOBIAN;
  1889.            continue;
  1890.            */
  1891.            fprintf(stderr, "ReadP3DFile: %s, Jacobian is unimplemented\n", p3dname);
  1892.            valid = 0;
  1893.            break;
  1894.          }
  1895.        if (!strcmp2(token, "NOJ") || !strcmp2(token, "NOJA") ||
  1896.            !strcmp2(token, "NOJAC") || !strcmp2(token, "NOJACO") ||
  1897.            !strcmp2(token, "NOJACOB") || !strcmp2(token, "NOJACOBI") ||
  1898.            !strcmp2(token, "NOJACOBIA") || !strcmp2(token, "NOJACOBIAN"))
  1899.          {
  1900.            /*
  1901.            if (isJacobian)
  1902.              {
  1903.                 fprintf(stderr, "ReadP3DFile: %s, double Jacobian spec.\n", p3dname);
  1904.                 valid = 0;
  1905.                 break; 
  1906.              }
  1907.            isJacobian = 1;
  1908.            check = P3D_NOJACOBIAN;
  1909.            */
  1910.            /* this is default */
  1911.            continue;
  1912.          }
  1913.  
  1914.        /* for blank spec */
  1915.        if (!strcmp2(token, "BL") || !strcmp2(token, "BLA") ||
  1916.            !strcmp2(token, "BLAN") || !strcmp2(token, "BLANK"))
  1917.          {
  1918.            if (isBlank && blank != P3D_BLANK)
  1919.              {
  1920.                 fprintf(stderr, "ReadP3DFile: %s, double Blank spec.\n", p3dname);
  1921.                 valid = 0;
  1922.                 break; /* get out of the while loop */ 
  1923.              }
  1924.            isBlank = 1;
  1925.            blank = P3D_BLANK;
  1926.            continue;
  1927.          }
  1928.        if (!strcmp2(token, "NOB") || !strcmp2(token, "NOBL") ||
  1929.            !strcmp2(token, "NOBLA") || !strcmp2(token, "NOBLAN") ||
  1930.            !strcmp2(token, "NOBLANK"))
  1931.          {
  1932.            if (isBlank && blank != P3D_NOBLANK)
  1933.              {
  1934.                 fprintf(stderr, "ReadP3DFile: %s, double Blank spec.\n", p3dname);
  1935.                 valid = 0;
  1936.                 break; /* get out of the while loop */ 
  1937.              }
  1938.            isBlank = 1;
  1939.            blank = P3D_NOBLANK;
  1940.            continue;
  1941.          }
  1942.  
  1943.        /* for Q file */
  1944.        if ( !strncmp2(token, "Q=", 2) )
  1945.          {
  1946.            p3d_set *curset = *sets + *nSets - 1; /* points to the current set */
  1947.            int nQ = curset->nQ; /* the number of Q files so far */
  1948.            p3d_q **qlistptr = &(curset->qlist); /* points to the qlist */
  1949.  
  1950.            /* if there is a filename */
  1951.            if (!token[2])
  1952.              {
  1953.                 fprintf(stderr, "ReadP3DFile: %s, missing Q file name.\n", p3dname);
  1954.                 valid = 0;
  1955.                 break; /* get out of the while loop */ 
  1956.              }
  1957.            
  1958.            /* alloc a storage for qlist and update the values */
  1959.            /* reset the array of qlist */
  1960.            if (nQ == 0)
  1961.              {
  1962.                *qlistptr = (p3d_q *) Alloc(sizeof(p3d_q));
  1963.              }
  1964.            else
  1965.              {
  1966.                *qlistptr = (p3d_q *) Realloc(*qlistptr, (nQ + 1) * sizeof(p3d_q));
  1967.              }
  1968.  
  1969.            /* if memory alloc is OK, then */
  1970.            if (!qlistptr)
  1971.              {
  1972.                fprintf(stderr, "unable to alloc memory\n");
  1973.                valid = 0;
  1974.                break;
  1975.              }
  1976.  
  1977.            /* use the old nQ value for array index */
  1978.            /* process the Q file spec */
  1979.            strcpy((*qlistptr)[nQ].name, token+2);
  1980.            (*qlistptr)[nQ].isMDataset = isMDataset;
  1981.            (*qlistptr)[nQ].fileType = fileType;
  1982.            (*qlistptr)[nQ].permutation = permutation; 
  1983.            (*qlistptr)[nQ].check = check;
  1984.            (*qlistptr)[nQ].jacobian = jacobian;
  1985.  
  1986.            /* increment the number of q file for the current p3d set */
  1987.            (curset->nQ)++;
  1988.  
  1989.            /* reset the flags exclusive to Q files */
  1990.            isFileType = 0;
  1991.            isPermutation = 0;
  1992.            isCheck = 0;
  1993.            isJacobian = 0;
  1994.  
  1995.            continue;
  1996.          }
  1997.  
  1998.        /* for FUNCTION file */
  1999.        if (!strncmp2(token, "F=", 2) || !strncmp2(token, "FU=", 3) ||
  2000.            !strncmp2(token, "FUN=", 4) || !strncmp2(token, "FUNC=", 5) ||
  2001.            !strncmp2(token, "FUNCT=", 6) || !strncmp2(token, "FUNCTI=", 7) ||
  2002.            !strncmp2(token, "FUNCTIO=", 8) || !strncmp2(token, "FUNCTION=", 9))
  2003.          {
  2004.            p3d_set *curset = *sets + *nSets - 1; /* points to the current set */
  2005.            int nFunc = curset->nFunc; /* the number of func files so far */
  2006.            p3d_func **funclistptr = &(curset->funclist); /* points to the funclist */
  2007.            int i = 0;
  2008.  
  2009.            /* search for the index of the file name */
  2010.            while(token[i++] != '=')
  2011.                ;
  2012.            /* if there is a filename */
  2013.            if (!token[i])
  2014.              {
  2015.                 fprintf(stderr, "ReadP3DFile: %s, missing FUNCTION file name.\n",
  2016.                          p3dname);
  2017.                 valid = 0;
  2018.                 break; /* get out of the while loop */ 
  2019.              }
  2020.            
  2021.            /* alloc a storage for qlist and update the values */
  2022.            /* reset the array of qlist */
  2023.            if (nFunc == 0)
  2024.              {
  2025.                *funclistptr = (p3d_func *) Alloc(sizeof(p3d_func));
  2026.              }
  2027.            else
  2028.              {
  2029.                *funclistptr = (p3d_func *) Realloc(*funclistptr,
  2030.                                       (nFunc + 1) * sizeof(p3d_func));
  2031.              }
  2032.  
  2033.            /* if memory alloc is OK, then */
  2034.            if (!funclistptr)
  2035.              {
  2036.                fprintf(stderr, "unable to alloc memory\n");
  2037.                valid = 0;
  2038.                break;
  2039.              }
  2040.  
  2041.            /* use the old nFunc value for array index */
  2042.            /* process the function file spec */
  2043.            strcpy((*funclistptr)[nFunc].name, token+i);
  2044.            (*funclistptr)[nFunc].fileType = fileType;
  2045.            (*funclistptr)[nFunc].permutation = permutation; 
  2046.  
  2047.            /* increment the number of function file for the current p3d set */
  2048.            (curset->nFunc)++;
  2049.  
  2050.            /* reset default values */
  2051.            isFileType = 0;
  2052.            isPermutation = 0;
  2053.  
  2054.            continue;
  2055.          }
  2056.  
  2057.  
  2058.        /* for GRID file */
  2059.        if (!strncmp2(token, "X=", 2) || !strncmp2(token, "XY=", 3) ||
  2060.            !strncmp2(token, "XYZ=", 4))
  2061.          {
  2062.            p3d_set *curset = *sets + *nSets - 1; /* points to the current set */
  2063.            p3d_grid **gridptr = &(curset->grid); /* points to the grid spec */
  2064.            int i = 0;
  2065.  
  2066.            while(token[i++] != '=')
  2067.                ;
  2068.  
  2069.            /* if there is a filename */
  2070.            if (!token[i])
  2071.              {
  2072.                fprintf(stderr, "ReadP3DFile: %s, missing XYZ file name.\n", p3dname);
  2073.                valid = 0;
  2074.                break; /* get out of the while loop */ 
  2075.              }
  2076.  
  2077.            /* it's an error to specify double XYZ for a single set */
  2078.            if (curset->isGrid)
  2079.              {
  2080.                fprintf(stderr, "ReadP3DFile: %s, double XYZ file spec.\n", p3dname);
  2081.                valid = 0;
  2082.                break; /* get out of the while loop */ 
  2083.              }
  2084.  
  2085.            /* alloc memory for grid file spec */
  2086.            *gridptr = (p3d_grid *) Alloc(sizeof(p3d_grid));
  2087.  
  2088.            /* if memory alloc is OK. then turn on the isGrid */
  2089.            if (!gridptr)
  2090.              {
  2091.                fprintf(stderr, "unable to alloc memory\n");
  2092.                valid = 0;
  2093.                break;
  2094.              }
  2095.  
  2096.            curset->isGrid = 1;
  2097.  
  2098.            /* process the grid file spec */
  2099.            strcpy((*gridptr)->name, token + i);
  2100.            (*gridptr)->fileType = fileType;
  2101.            (*gridptr)->permutation = permutation; 
  2102.            (*gridptr)->blank = blank;
  2103.  
  2104.            /* reset default values */
  2105.            isFileType = 0;
  2106.            isPermutation = 0;
  2107.            isBlank = 0;
  2108.  
  2109.            continue;
  2110.          }
  2111.  
  2112.        if (!strcmp2(token, "READ"))
  2113.          {
  2114.            p3d_set *curset = *sets + *nSets - 1;
  2115.            p3d_set *temp;
  2116.            int k;
  2117.  
  2118.            /* check the set processed so far */
  2119.            if (!curset->isGrid)
  2120.              {
  2121.                fprintf(stderr, "ReadP3DFile: in %s XYZ file missing.\n", p3dname);
  2122.                valid = 0;
  2123.                break;
  2124.              }
  2125.            else if(!curset->nQ && !curset->nFunc)
  2126.              {
  2127.                fprintf(stderr, "ReadP3DFile: %s, no Q or function file.\n", p3dname);
  2128.                valid = 0;
  2129.                break;
  2130.              }
  2131.  
  2132.            curset->ndim = ndim;
  2133.            curset->isMGrid = isMGrid;
  2134.  
  2135.            /* another set begins here */
  2136.            if (!(temp = (p3d_set *) Alloc(sizeof(p3d_set) * (*nSets + 1))))
  2137.              {
  2138.                fprintf(stderr, "Unable to alloc memory, aborted\n");
  2139.                valid = 0;
  2140.                break; /* get out of the while loop */ 
  2141.              }
  2142.            for (k = 0; k < *nSets; k++)
  2143.              {
  2144.                temp[k] = (*sets)[k];
  2145.              }
  2146.            Free(*sets); 
  2147.            *sets = temp;
  2148.            (*nSets)++;
  2149.            curset = *sets + *nSets - 1;
  2150.  
  2151.            /* initialize some values */
  2152.            curset->ndim = 3;
  2153.            curset->isMGrid = 0;
  2154.            curset->isGrid = 0;
  2155.            curset->grid = NULL;
  2156.            curset->nQ = 0;
  2157.            curset->qlist = NULL;
  2158.            curset->nFunc = 0;
  2159.            curset->funclist = NULL;
  2160.  
  2161.            /* reset default values */
  2162.            isNdim = 0, ndim = 3;
  2163.            isMDataset = 0;
  2164.            isMGrid = 0;
  2165.            isFileType = 0, fileType = P3D_DAT;
  2166.            isPermutation = 0, permutation = P3D_WHOLE;
  2167.            isCheck = 0, check = P3D_NOCHECK;
  2168.            isJacobian = 0, jacobian = P3D_NOJACOBIAN;
  2169.            isBlank = 0, blank = P3D_NOBLANK;
  2170.  
  2171.            continue;
  2172.          }
  2173.        
  2174.        fprintf(stderr, "ReadP3DFile: %s, unknown spec %s.\n", p3dname, token);
  2175.        valid = 0;
  2176.        break; /* get out of the while loop */ 
  2177.     } /* while */
  2178.  
  2179.    fclose(inFile);
  2180.  
  2181.    if (!valid || status == -1)
  2182.      {
  2183.        fprintf(stderr, "ReadP3DFile: %s, file format error.\n", p3dname);
  2184.        return 0;
  2185.      }
  2186.  
  2187.    /* also, check current set, about nQ, nFunc, and isGrid, to see
  2188.       if there is any error */
  2189.    if ((*sets + *nSets - 1)->isGrid == 0)
  2190.      {
  2191.        fprintf(stderr, "ReadP3DFile: %s, XYZ missing in number %d read command.\n", p3dname, *nSets);
  2192.        return 0;
  2193.      }
  2194.    else if((*sets + *nSets - 1)->nQ == 0 && (*sets + *nSets - 1)->nFunc == 0)
  2195.      {
  2196.        fprintf(stderr, "ReadP3DFile: %s, no Q or function file.\n", p3dname);
  2197.        return 0;
  2198.      }
  2199.    else
  2200.      {
  2201.        (*sets + *nSets - 1)->ndim = ndim;
  2202.        (*sets + *nSets - 1)->isMGrid = isMGrid;
  2203.      }
  2204.  
  2205.    return 1;
  2206. } /* parseP3DMeta */
  2207.  
  2208. /* using row major order arrangement */
  2209. #ifdef PROTO
  2210. static void putIblank(int iblank, int **ipp, long *indices, long *dims,
  2211.         int ndim, int igrid) 
  2212. #else
  2213. static void putIblank(iblank, ipp, indices, dims, ndim, igrid) 
  2214. int iblank;
  2215. int **ipp;
  2216. long *indices;
  2217. long *dims;
  2218. int ndim;  /* rank */
  2219. int igrid;
  2220. #endif
  2221. /* return 0 when failed */
  2222. {
  2223.   int i;
  2224.   long offset;
  2225.  
  2226.   offset = indices[ndim - 1];
  2227.   for (i = ndim - 1; i > 0; i--)
  2228.      {
  2229.        offset = offset * dims[i - 1] + indices[i - 1];
  2230.      }
  2231.   ipp[igrid][offset] = iblank;
  2232.  
  2233.   return;
  2234. } /* putIblank */
  2235.  
  2236. #ifdef PROTO
  2237. static void getIblank(int *iblank, int **ipp, long *indices,
  2238.                                   long *dims, int ndim, int igrid) 
  2239. #else
  2240. static void getIblank(iblank, ipp, indices, dims, ndim, igrid)
  2241. int *iblank;
  2242. int **ipp;
  2243. long *indices;
  2244. long *dims;
  2245. int ndim;
  2246. int igrid;
  2247. #endif
  2248. {
  2249.   int i;
  2250.   long offset;
  2251.  
  2252.   offset = indices[ndim - 1];
  2253.   for (i = ndim - 1; i > 0; i--)
  2254.      {
  2255.        offset = offset * dims[i - 1] + indices[i - 1];
  2256.      }
  2257.   *iblank = ipp[igrid][offset];
  2258.  
  2259.   return;
  2260. } /* getIblank */
  2261.  
  2262. #ifdef PROTO
  2263. static void freeiblank(int ngrid, int **ipp)
  2264. #else
  2265. static void freeiblank(ngrid, ipp)
  2266. int ngrid;
  2267. int **ipp;
  2268. #endif
  2269. {
  2270.   int i;
  2271.  
  2272.   for (i = 0; i < ngrid; i++)
  2273.      Free(ipp[i]);
  2274.   Free(ipp);
  2275.  
  2276.   return;
  2277. }
  2278.  
  2279. #ifdef PROTO
  2280. static void freedims(int ngrid, long **dimsptr)
  2281. #else
  2282. static void freedims(ngrid, dimsptr)
  2283. int ngrid;
  2284. long **dimsptr;
  2285. #endif
  2286. {
  2287.   int i;
  2288.  
  2289.   for (i = 0; i < ngrid; i++)
  2290.      Free(dimsptr[i]);
  2291.   Free(dimsptr);
  2292.  
  2293.   return;
  2294. }
  2295.  
  2296. /*********************formatted p3d reading routines, start *********************/
  2297.  
  2298. #ifdef PROTO
  2299. static int RdFmtXYZ(int ndim, int isMGrid, long ***dimsptr,
  2300.                                p3d_grid *gridptr, ObjPtr **formptr, int ***ipp)
  2301. #else
  2302. static int RdFmtXYZ(ndim, isMGrid, dimsptr, gridptr, formptr, ipp)
  2303. int ndim;
  2304. int isMGrid;
  2305. long ***dimsptr;
  2306. p3d_grid *gridptr;
  2307. ObjPtr **formptr;
  2308. int ***ipp;
  2309. #endif
  2310. {
  2311.   FILE *inFile;
  2312.   int ngrid;     /* number of grids */
  2313.   long *indices;
  2314.   real *bounds;
  2315.   int i, j;
  2316.   int n;
  2317.   char *name;  /* dataset name starting point */
  2318.   int namelen;
  2319.   int repeat;
  2320.   double dval;
  2321.  
  2322.   if (!(inFile = fopen(gridptr->name, "r")))
  2323.     {
  2324.       fprintf(stderr, "Unable to open %s, aborted\n", gridptr->name);
  2325.       return 0;
  2326.     }
  2327.  
  2328.   /* determine the dataset name now */
  2329.   for (i = (int) strlen(gridptr->name) - 1; i >= 0 && gridptr->name[i] != '/'; i--)
  2330.      ;
  2331.   if (i == -1 || gridptr->name[i] == '/')
  2332.     i++;
  2333.   name = gridptr->name + i;
  2334.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  2335.      ;
  2336.  
  2337.   if (isMGrid)
  2338.     {
  2339.       SkipBlanksAndCommas(inFile);
  2340.       if (fscanf(inFile, "%d", &ngrid) != 1 || ngrid <= 0)
  2341.         {
  2342.           fprintf(stderr, "%s format error, invalid NGRID value\n", gridptr->name);
  2343.           fclose(inFile);
  2344.           return 0;
  2345.         } 
  2346.     }
  2347.   else
  2348.     {
  2349.       ngrid = 1;
  2350.     }
  2351.   
  2352.   /* alloc dimsptr and indices storage */
  2353.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  2354.     {
  2355.       fprintf(stderr, "Unable to alloc memory\n");
  2356.       fclose(inFile);
  2357.       return 0;
  2358.     }
  2359.   if (!(*dimsptr = (long **) Alloc(ngrid * sizeof(long *))))
  2360.     {
  2361.       fprintf(stderr, "Unable to alloc memory\n");
  2362.       Free(indices);
  2363.       fclose(inFile);
  2364.       return 0;
  2365.     }
  2366.   repeat = 0;
  2367.   for (i = 0; i < ngrid; i++)
  2368.     {
  2369.       if (!((*dimsptr)[i] = (long *) Alloc(ndim * sizeof(long))))
  2370.         {
  2371.           fprintf(stderr, "Unable to alloc memory\n");
  2372.           Free(indices);
  2373.           freedims(i, *dimsptr);
  2374.           fclose(inFile);
  2375.           return 0;
  2376.         }
  2377.  
  2378.       for (j = 0; j < ndim; j++)
  2379.          {
  2380.            if (!repeat)
  2381.              {
  2382.                 if ((repeat = getValue(inFile, &dval)) <= 0 || (int) dval <= 0)
  2383.                   {
  2384.                     fprintf(stderr, "Invalid dimension value in %s\n",
  2385.                             gridptr->name);
  2386.                     Free(indices);
  2387.                     freedims(i, *dimsptr);
  2388.                     fclose(inFile);
  2389.                     return 0;
  2390.                   }
  2391.              }
  2392.            repeat--;
  2393.            *((*dimsptr)[i] + j) = (long) dval; 
  2394.          }
  2395.     }
  2396.  
  2397.   /* alloc strorage for pointer to bounds */
  2398.   if (!(bounds = (real *) Alloc(ndim * 2 * sizeof(real))))
  2399.     {
  2400.       fprintf(stderr, "Unable to alloc memory\n");
  2401.       Free(indices);
  2402.       freedims(ngrid, *dimsptr);
  2403.       fclose(inFile);
  2404.       return 0;
  2405.     }
  2406.   /* alloc storage for pointers to dataforms */
  2407.   if (!(*formptr = (ObjPtr *) Alloc(ngrid * sizeof(ObjPtr))))
  2408.     {
  2409.       fprintf(stderr, "Unable to alloc memory\n");
  2410.       Free(indices);
  2411.       Free(bounds);
  2412.       freedims(ngrid, *dimsptr);
  2413.       fclose(inFile);
  2414.       return 0;
  2415.     }
  2416.  
  2417.   /* alloc iblankptr if IBLANK */
  2418.   if (gridptr->blank == P3D_BLANK)
  2419.     {
  2420.       if (!(*ipp = (int **) Alloc(ngrid * sizeof(int *))))
  2421.         {
  2422.           fprintf(stderr, "Unable to alloc memory\n");
  2423.           Free(indices);
  2424.           Free(bounds);
  2425.           Free(*formptr);
  2426.           freedims(ngrid, *dimsptr);
  2427.           fclose(inFile);
  2428.           return 0;
  2429.         }
  2430.  
  2431.     }
  2432.  
  2433.   repeat = 0;
  2434.   for (n = 0; n < ngrid; n++)
  2435.     {
  2436.       long size, count; /* number of grid values promised */
  2437.       ObjPtr dataset;
  2438.       char datasetname[MAXNAMELEN];
  2439.       int whichcomponent;
  2440.  
  2441.       strncpy(datasetname, name, namelen);
  2442.       datasetname[namelen] = '\0';
  2443.       if (ngrid > 1)
  2444.         {
  2445.           sprintf(datasetname + namelen, "%d", n + 1);
  2446.         }
  2447.  
  2448.       /* setup a dataset for the grid data */
  2449.       dataset = NewStructuredDataset(datasetname, ndim, (*dimsptr)[n], ndim);
  2450.       if (!SetCurField(FIELD1, dataset))
  2451.         {
  2452.           fprintf(stderr, "Unable to SetCurField\n");
  2453.           Free(indices);
  2454.           Free(bounds);
  2455.           Free(*formptr);
  2456.           freedims(ngrid, *dimsptr);
  2457.           if (gridptr->blank == P3D_BLANK)
  2458.             {
  2459.               freeiblank(n-1, *ipp);
  2460.             }
  2461.           fclose(inFile);
  2462.           return 0;
  2463.         }
  2464.  
  2465.       size = 1;
  2466.       for (i = 0; i < ndim; i++)
  2467.          {
  2468.            size *= (*dimsptr)[n][i];
  2469.          }
  2470.       size *= (long)(ndim + (gridptr->blank == P3D_BLANK));
  2471.  
  2472.       if (gridptr->blank == P3D_BLANK)
  2473.         {
  2474.           /* alloc storage for the iblank arrays */
  2475.           if (!((*ipp)[n] = (int *) Alloc(size * sizeof(int))))
  2476.             {
  2477.               fprintf(stderr, "Unable to alloc memory\n");
  2478.               freeiblank(n, *ipp);
  2479.               Free(indices);
  2480.               Free(bounds);
  2481.               Free(*formptr);
  2482.               freedims(ngrid, *dimsptr);
  2483.               fclose(inFile);
  2484.               return 0;
  2485.             }
  2486.         }
  2487.  
  2488.       for (i = 0; i < ndim; i++)
  2489.          {
  2490.            indices[i] = 0;
  2491.          }
  2492.       whichcomponent = 0;
  2493.       /* read in all the grid values */
  2494.       for (count = 0; count < size; count++)
  2495.         {
  2496.           int first = 1;
  2497.  
  2498.           if (!repeat)
  2499.             {
  2500.               if ((repeat = getValue(inFile, &dval)) <= 0)
  2501.                 {
  2502.                   fprintf(stderr, "Error reading %s, invalid numeric value\n",
  2503.                            gridptr->name);
  2504.                   Free(indices);
  2505.                   Free(bounds);
  2506.                   Free(*formptr);
  2507.                   freedims(ngrid, *dimsptr);
  2508.                   if (gridptr->blank == P3D_BLANK)
  2509.                     {
  2510.                       freeiblank(n + 1, *ipp);
  2511.                     }
  2512.                   fclose(inFile);
  2513.                   return 0;
  2514.                 }
  2515.             }
  2516.  
  2517.           repeat--;
  2518.           if (whichcomponent < ndim) /* read grid values */
  2519.             {
  2520.               PutFieldComponent(FIELD1, whichcomponent, indices, (real) dval);  
  2521.               for (j = 0; first && j < ndim; j++)
  2522.                  first = indices[j] == 0;
  2523.  
  2524.               if (first)
  2525.                 {
  2526.                   bounds[whichcomponent * 2] = (real) dval;
  2527.                   bounds[whichcomponent * 2 + 1] = (real) dval;
  2528.                 }
  2529.               else
  2530.                 {
  2531.                   if (dval < bounds[whichcomponent * 2])
  2532.                     bounds[whichcomponent * 2] = (real) dval;
  2533.                   else if (dval > bounds[whichcomponent * 2 + 1])
  2534.                     bounds[whichcomponent * 2 + 1] = (real) dval;
  2535.                 }
  2536.             }
  2537.           else /* read iblank values */
  2538.             {
  2539.               putIblank((int) dval, *ipp, indices, (*dimsptr)[n], ndim, n);
  2540.             }
  2541.           /* adjust the indices and whichcomponent */
  2542.           if (gridptr->permutation == P3D_WHOLE)
  2543.             {
  2544.               for (j = 0; j < ndim; j++)
  2545.                  {
  2546.                    indices[j]++;
  2547.                    if (indices[j] == (*dimsptr)[n][j])
  2548.                      {
  2549.                        indices[j] = 0;
  2550.                        continue;
  2551.                      }
  2552.                    break;
  2553.                  }
  2554.               if (j == ndim) /* all indices are set to zero now */ 
  2555.                 {
  2556.                   whichcomponent++;
  2557.                 }
  2558.             }
  2559.           else  /* P3D_PLANE */
  2560.             {
  2561.               for (j = 0; j < ndim - 1; j++)
  2562.                  {
  2563.                    indices[j]++;
  2564.                    if (indices[j] == (*dimsptr)[n][j])
  2565.                      {
  2566.                        indices[j] = 0;
  2567.                        continue;
  2568.                      }
  2569.                    break;
  2570.                  }
  2571.               if (j == ndim - 1)  
  2572.                 {
  2573.                   whichcomponent++;
  2574.                   if ((gridptr->blank == P3D_BLANK) ? 
  2575.                      (whichcomponent == ndim + 1) : (whichcomponent == ndim))
  2576.                   /* leave a space for reading iblanks */
  2577.                     {
  2578.                       whichcomponent = 0;
  2579.                       indices[j]++;
  2580.                     }
  2581.                 }
  2582.             }
  2583.         } /* for read grid values */
  2584.  
  2585.      (*formptr)[n] = NewCurvilinearDataForm(datasetname, ndim, (*dimsptr)[n],
  2586.                     bounds, dataset);
  2587.  
  2588.     } /* outer for */
  2589.  
  2590.   /* indices and bounds are no longer of any use */
  2591.   Free(indices);
  2592.   Free(bounds);
  2593.   fclose(inFile);
  2594.   return n;
  2595. } /* RdFmtXYZ */
  2596.  
  2597. #ifdef PROTO
  2598. static int RdFmtSol(int ndim, int isMGrid, int ngrid, int blank,
  2599.                  long **dimsptr, p3d_q *qptr, ObjPtr *formptr, int **ipp)
  2600. #else
  2601. static int RdFmtSol(ndim, isMGrid, ngrid, blank, dimsptr, qptr, formptr, ipp)
  2602. int ndim;
  2603. int isMGrid;
  2604. int ngrid;
  2605. int blank;
  2606. long **dimsptr;
  2607. p3d_q *qptr;
  2608. ObjPtr *formptr;
  2609. int **ipp;
  2610. #endif
  2611. {
  2612.   FILE *inFile;
  2613.   int i, j, n;
  2614.   char *name;
  2615.   long *indices;
  2616.   int numsets; /* number of datasets read so far, for multiple datasets in a file*/
  2617.   int namelen; /* the length of dataset name going to be used */
  2618.   int repeat;  /* to deal with <num>*<value> format */
  2619.   double dval;
  2620.  
  2621.   if (!(inFile = fopen(qptr->name, "r")))
  2622.     {
  2623.       fprintf(stderr, "Unable to open %s, aborted\n", qptr->name);
  2624.       return 0;
  2625.     }
  2626.  
  2627.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  2628.     {
  2629.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  2630.       fclose(inFile);
  2631.       return 0;
  2632.     }
  2633.  
  2634.   /* determine the dataset name from the file name */
  2635.   for (i = (int) strlen(qptr->name) - 1; i >= 0 && qptr->name[i] != '/'; i--)
  2636.      ;
  2637.   if (i == -1 || qptr->name[i] == '/')
  2638.     i++;
  2639.   name = qptr->name + i;
  2640.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  2641.      ;
  2642.  
  2643.   numsets = 0;
  2644.  
  2645.   do
  2646.     {
  2647.       repeat = 0;
  2648.       /* get ngrid value if it is a MGRID q file */
  2649.       if (isMGrid)
  2650.         {
  2651.            repeat = getValue(inFile, &dval);
  2652.            if (repeat == EOF)
  2653.              {
  2654.               if (numsets == 0)
  2655.                  {
  2656.                   fprintf(stderr, "Error in %s: empty file\n", qptr->name);
  2657.                  }
  2658.               fclose(inFile);
  2659.               Free(indices);
  2660.               return (numsets ? ngrid : 0);
  2661.              }
  2662.           else if (repeat == 0 || (n = (int) dval) <= 0)
  2663.              {
  2664.               fprintf(stderr, "Error in %s: Invalid NGRID spec.\n", qptr->name);
  2665.               fclose(inFile);
  2666.               Free(indices);
  2667.               return 0;
  2668.              }
  2669.           repeat--;
  2670.         }
  2671.       else
  2672.         {
  2673.           n = 1;
  2674.         }
  2675.  
  2676.       /* if ngrid != this file's ngrid then error */
  2677.       if (ngrid != n)
  2678.         {
  2679.           fprintf(stderr, "Error in %s: NGRID value does not match its grid\n",
  2680.                                                              qptr->name);
  2681.           fclose(inFile);
  2682.           Free(indices);
  2683.           return 0;
  2684.           
  2685.         }
  2686.       /* if dims does not match the dims read in XYZ file then error */
  2687.       for (n = 0; n < ngrid; n++)
  2688.         {
  2689.           for (i = 0; i < ndim; i++)
  2690.             {
  2691.               if (!repeat)
  2692.                 {
  2693.                   if ((repeat = getValue(inFile, &dval)) <= 0)
  2694.                     {
  2695.                       /* if this is not multiple grid, then the first
  2696.                       dimension value is used to indicate whether there is
  2697.                       more datasets or not */
  2698.                       if (!isMGrid && i == 0 && repeat == EOF)
  2699.                         {
  2700.                          if (numsets == 0)
  2701.                            {
  2702.                             fprintf(stderr, "Error in %s: empty file\n", qptr->name);
  2703.                            }
  2704.                           fclose(inFile);
  2705.                           Free(indices);
  2706.                           return (numsets ? ngrid : 0);
  2707.                         }
  2708.                       fprintf(stderr, "Invalid dimension value in %s\n", qptr->name);
  2709.                       fclose(inFile);
  2710.                       Free(indices);
  2711.                       return 0;
  2712.                     }
  2713.                 }
  2714.               repeat--;
  2715.               if ((long) dval != dimsptr[n][i])
  2716.                 {
  2717.                   fprintf(stderr,
  2718.                   "Dimensions mismatch between %s and its XYZ file\n", qptr->name);
  2719.                   fclose(inFile);
  2720.                   Free(indices);
  2721.                   return 0;
  2722.                 }
  2723.             }
  2724.         }
  2725.     
  2726.       for (n = 0; n < ngrid; n++)
  2727.         {
  2728.           long size, count; /* number of data values promised */
  2729.           ObjPtr density, pressure, dataset;
  2730.           char datasetname[MAXNAMELEN];
  2731.           char densityname[MAXNAMELEN];
  2732.           char pressurename[MAXNAMELEN];
  2733.           double time;
  2734.           int whichcomponent;
  2735.           ObjPtr whichdataset;
  2736.     
  2737.           /* process the FAMACH..TIME line */
  2738.           for (i = 1; i <= 4; i++)
  2739.             {
  2740.               if (!repeat)
  2741.                 {
  2742.                   if ((repeat = getValue(inFile, &dval)) <= 0)
  2743.                     {
  2744.                      fprintf(stderr,
  2745.                      "Error in %s: invalid FSMACH,ALPHA,RE,TIME line\n", qptr->name);
  2746.                      Free(indices);
  2747.                      fclose(inFile);
  2748.                      return 0;
  2749.                     }
  2750.                 }
  2751.               repeat--;
  2752.               if (i == 4)
  2753.                 {
  2754.                   time = dval;
  2755.                 }
  2756.             }
  2757.     
  2758.           /* setup datasets */
  2759.           strncpy(datasetname, name, namelen);
  2760.           datasetname[namelen] = '\0';
  2761.           if (ngrid > 1)
  2762.             {
  2763.           /* the ordering starts from 1 */
  2764.               sprintf(datasetname + namelen, "%d", n + 1);
  2765.             }
  2766.           sprintf(densityname, "%s%s", datasetname, " Density");
  2767.           sprintf(pressurename, "%s%s", datasetname, " Pressure");
  2768.           density = NewStructuredDataset(densityname, ndim, dimsptr[n], 0);
  2769.           pressure = NewStructuredDataset(pressurename, ndim, dimsptr[n], 0);
  2770.           dataset = NewStructuredDataset(datasetname, ndim, dimsptr[n], ndim);
  2771.  
  2772.           size = 1;
  2773.           for (i = 0; i < ndim; i++)
  2774.              size *= dimsptr[n][i];
  2775.           size *= ndim + 2;  /* includes density and pressure */
  2776.           for (i = 0; i < ndim; i++)
  2777.              {
  2778.                indices[i] = 0;
  2779.              }
  2780.           whichcomponent = 0;
  2781.           whichdataset = density;
  2782.     
  2783.           /* read in all the q values */
  2784.           for (count = 0; count < size; count++)
  2785.             {
  2786.               int iblank = 1;
  2787.     
  2788.               if (!repeat)
  2789.                 {
  2790.                   if ((repeat = getValue(inFile, &dval)) <= 0)
  2791.                     {
  2792.                       fprintf(stderr, "Error reading %s, invalid numeric value\n",
  2793.                                   qptr->name);
  2794.                       Free(indices);
  2795.                       fclose(inFile);
  2796.                       return 0;
  2797.                     }
  2798.                 }
  2799.               repeat--;
  2800.               /* do check of density and pressure, and iblanks here */
  2801.               if (qptr->check == P3D_CHECK && (whichdataset == pressure ||
  2802.                   whichdataset == density))
  2803.                 {
  2804.                   if (dval < 0.0)
  2805.             {
  2806.               fprintf(stderr, "Negative %s value %f found in Q file %s and is converted to %f\n", whichdataset == pressure ? "pressure" : "density", dval, qptr->name, 0.0);
  2807.                       dval = 0.0;
  2808.                     }
  2809.                 }
  2810.               if (blank == P3D_BLANK)
  2811.                 {
  2812.                    getIblank(&iblank, ipp, indices, dimsptr[n], ndim, n);
  2813.                 }
  2814.     
  2815.               if (!SetCurField(FIELD1, whichdataset))
  2816.                 {
  2817.                   fprintf(stderr, "Unable to SetCurField\n");
  2818.                   Free(indices);
  2819.                   fclose(inFile);
  2820.                   return 0;
  2821.                 }
  2822.     
  2823.               PutFieldComponent(FIELD1, whichcomponent, indices, 
  2824.                           iblank == 0 ? (real) missingData : (real) dval);  
  2825.     
  2826.               /* adjust the indices, whichdataset, and whichcomponent */
  2827.               if (qptr->permutation == P3D_WHOLE)
  2828.                 {
  2829.                   for (j = 0; j < ndim; j++)
  2830.                      {
  2831.                        indices[j]++;
  2832.                        if (indices[j] == dimsptr[n][j])
  2833.                          {
  2834.                            indices[j] = 0;
  2835.                            continue;
  2836.                          }
  2837.                        break;
  2838.                      }
  2839.                   if (j == ndim) /* all indices are set to zero now */ 
  2840.                     {
  2841.                       if (whichdataset == dataset)
  2842.                         {
  2843.                           whichcomponent++;
  2844.                         }
  2845.                       else
  2846.                         {
  2847.                           if (whichdataset == density)
  2848.                             whichdataset = pressure;
  2849.                           else
  2850.                             whichdataset = dataset;
  2851.                         }
  2852.                     }
  2853.                 }
  2854.               else  /* P3D_PLANE */
  2855.                 {
  2856.                   for (j = 0; j < ndim - 1; j++)
  2857.                      {
  2858.                        indices[j]++;
  2859.                        if (indices[j] == dimsptr[n][j])
  2860.                          {
  2861.                            indices[j] = 0;
  2862.                            continue;
  2863.                          }
  2864.                        break;
  2865.                      }
  2866.                   if (j == ndim - 1)  
  2867.                     {
  2868.                       if (whichdataset == density)
  2869.                         whichdataset = pressure;
  2870.                       else if (whichdataset == pressure)
  2871.                         whichdataset = dataset;
  2872.                       else /* whichdataset == dataset */
  2873.                         {
  2874.                           whichcomponent++;
  2875.                           if (whichcomponent == ndim)
  2876.                             {
  2877.                               whichdataset = density;
  2878.                               whichcomponent = 0;
  2879.                               indices[j]++;
  2880.                             }
  2881.                         }
  2882.                     }
  2883.                 }
  2884.             } 
  2885.     
  2886.           SetDatasetForm(density, formptr[n]);
  2887.           SetDatasetForm(pressure, formptr[n]);
  2888.           SetDatasetForm(dataset, formptr[n]);
  2889.           if (qptr->isMDataset)
  2890.             {
  2891.               RegisterTimeDataset(density, (real) time);
  2892.               RegisterTimeDataset(pressure, (real) time);
  2893.               RegisterTimeDataset(dataset, (real) time);
  2894.             }
  2895.           else
  2896.             {
  2897.               RegisterDataset(density);
  2898.               RegisterDataset(pressure);
  2899.               RegisterDataset(dataset);
  2900.             }
  2901.         } /* outer for */
  2902.       numsets++;
  2903.     } while(qptr->isMDataset);
  2904.  
  2905.   Free(indices);
  2906.   fclose(inFile);
  2907.   return ngrid;
  2908. } /* RdFmtSol */
  2909.  
  2910. #ifdef PROTO
  2911. static int RdFmtFun(int ndim, int isMGrid, int ngrid, int blank,
  2912.                 long **dimsptr, p3d_func *funcptr, ObjPtr *formptr, int **ipp)
  2913. #else
  2914. static int RdFmtFun(ndim, isMGrid, ngrid, blank, 
  2915.                               dimsptr, funcptr, formptr, ipp)
  2916. int ndim;
  2917. int isMGrid;
  2918. int ngrid;
  2919. int blank;
  2920. long **dimsptr;
  2921. p3d_func *funcptr;
  2922. ObjPtr *formptr;
  2923. int **ipp;
  2924. #endif
  2925. {
  2926.   FILE *inFile;/* the current function file being processed */
  2927.   int i, j, n;
  2928.   char *name; /* points to the portion of function file name after the last / */
  2929.   long *indices;
  2930.   int *nvar;  /* number of variables of this function */
  2931.   int namelen;
  2932.   int repeat;
  2933.   double dval;
  2934.  
  2935.   if (!(inFile = fopen(funcptr->name, "r")))
  2936.     {
  2937.       fprintf(stderr, "Unable to open %s, aborted\n", funcptr->name);
  2938.       return 0;
  2939.     }
  2940.   if (isMGrid)
  2941.     {
  2942.       SkipBlanksAndCommas(inFile);
  2943.       if (fscanf(inFile, "%d", &n) != 1 || n <= 0)
  2944.         {
  2945.           fprintf(stderr, "Error in %s: invalid NGRID value\n", funcptr->name);
  2946.           fclose(inFile);
  2947.           return 0;
  2948.         }
  2949.     }
  2950.   else
  2951.     {
  2952.       n = 1;
  2953.     }
  2954.   /* if ngrid != this file's ngrid then error */
  2955.   if (ngrid != n)
  2956.     {
  2957.       fprintf(stderr, "Error in %s: NGRID value does not match its grid\n",
  2958.                                                                  funcptr->name);
  2959.       fclose(inFile);
  2960.       return 0;
  2961.     }
  2962.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  2963.     {
  2964.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  2965.       fclose(inFile);
  2966.       return 0;
  2967.     }
  2968.   if (!(nvar = (int *) Alloc(ngrid * sizeof(int))))
  2969.     {
  2970.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  2971.       Free(indices);
  2972.       fclose(inFile);
  2973.       return 0;
  2974.     }
  2975.  
  2976.   /* if dims does not match the dims read in XYZ file then error */
  2977.   repeat = 0;
  2978.   for (n = 0; n < ngrid; n++)
  2979.     {
  2980.       for (i = 0; i < ndim + 1; i++) /* nvar follows dimension spec's */
  2981.         {
  2982.           if (!repeat)
  2983.             {
  2984.               if ((repeat = getValue(inFile, &dval)) <= 0)
  2985.                 {
  2986.                   fprintf(stderr, "Invalid dim/nvar value in %s\n", funcptr->name);
  2987.                   Free(indices);
  2988.                   Free(nvar);
  2989.                   fclose(inFile);
  2990.                   return 0;
  2991.                 }
  2992.             }
  2993.           repeat--;
  2994.           if (i < ndim)  /* the value is dimension */
  2995.             {
  2996.               if ((long) dval != dimsptr[n][i])
  2997.                 {
  2998.                   fprintf(stderr, "Dimensions mismatch between %s"
  2999.                                       " and its XYZ file\n", funcptr->name);
  3000.                   Free(indices);
  3001.                   Free(nvar);
  3002.                   fclose(inFile);
  3003.                   return 0;
  3004.                 }
  3005.             }
  3006.           else /* i = ndim, the value is of nvar */
  3007.             {
  3008.               nvar[n] = (int) dval;
  3009.             }
  3010.         }
  3011.     }
  3012.  
  3013.   /* determine the dataset name */
  3014.   for (i = (int) strlen(funcptr->name) - 1; i >= 0 && funcptr->name[i] != '/'; i--)
  3015.      ;
  3016.   if (i == -1 || funcptr->name[i] == '/')
  3017.     i++;
  3018.   name = funcptr->name + i;
  3019.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  3020.      ;
  3021.  
  3022.   repeat = 0;
  3023.   for (n = 0; n < ngrid; n++)
  3024.     {
  3025.       long size, count; /* number of data values promised */
  3026.       ObjPtr dataset;
  3027.       char datasetname[MAXNAMELEN];
  3028.       int whichcomponent;
  3029.  
  3030.       strncpy(datasetname, name, namelen);
  3031.       datasetname[namelen] = '\0';
  3032.       /* number the datasets only if its number of grids is than 1 */
  3033.       if (ngrid > 1)
  3034.         {
  3035.       /* ordering starts from 1 */
  3036.           sprintf(datasetname + namelen, "%d", n + 1);
  3037.         }
  3038.  
  3039.       size = 1;
  3040.       for (i = 0; i < ndim; i++)
  3041.          size *= dimsptr[n][i];
  3042.       size *= nvar[n];
  3043.       for (i = 0; i < ndim; i++)
  3044.          {
  3045.            indices[i] = 0;
  3046.          }
  3047.  
  3048.       /* setup datasets */
  3049.       dataset = NewStructuredDataset(datasetname, ndim, dimsptr[n],
  3050.                                               nvar[n] == 1 ? 0 : nvar[n]);
  3051.       if (!SetCurField(FIELD1, dataset))
  3052.         {
  3053.           fprintf(stderr, "Unable to SetCurField\n");
  3054.           Free(indices);
  3055.           Free(nvar);
  3056.           fclose(inFile);
  3057.           return 0;
  3058.         }
  3059.       whichcomponent = 0;
  3060.  
  3061.       /* read in all the func values */
  3062.       for (count = 0; count < size; count++)
  3063.         {
  3064.           int iblank = 1;
  3065.      
  3066.           if (!repeat)
  3067.             {
  3068.               if ((repeat = getValue(inFile, &dval)) <= 0)
  3069.                 {
  3070.                   fprintf(stderr, "Error reading %s, invalid numeric value\n",
  3071.                                 funcptr->name);
  3072.                   Free(indices);
  3073.                   Free(nvar);
  3074.                   fclose(inFile);
  3075.                   return 0;
  3076.                 }
  3077.             }
  3078.           repeat--;
  3079.           if (blank == P3D_BLANK)
  3080.             {
  3081.                getIblank(&iblank, ipp, indices, dimsptr[n], ndim, n);
  3082.             }
  3083.     
  3084.           PutFieldComponent(FIELD1, whichcomponent, indices, 
  3085.                   iblank == 0 ? (real) missingData : (real) dval);  
  3086.  
  3087.           /* adjust the indices and whichcomponent */
  3088.           if (funcptr->permutation == P3D_WHOLE)
  3089.             {
  3090.               for (j = 0; j < ndim; j++)
  3091.                 {
  3092.                   indices[j]++;
  3093.                   if (indices[j] == dimsptr[n][j])
  3094.                     {
  3095.                       indices[j] = 0;
  3096.                       continue;
  3097.                     }
  3098.                   break;
  3099.                 }
  3100.               if (j == ndim) /* all indices are set to zero now */
  3101.                 {
  3102.                   whichcomponent++;
  3103.                 }
  3104.             }
  3105.           else  /* P3D_PLANE */
  3106.             {
  3107.               for (j = 0; j < ndim - 1; j++)
  3108.                 {
  3109.                   indices[j]++;
  3110.                   if (indices[j] == dimsptr[n][j])
  3111.                     {
  3112.                       indices[j] = 0;
  3113.                       continue;
  3114.                     }
  3115.                   break;
  3116.                 }
  3117.               if (j == ndim - 1 && ++whichcomponent == nvar[n])
  3118.                 {
  3119.                   indices[j]++;
  3120.                   whichcomponent = 0;
  3121.                 }
  3122.             }
  3123.         }     /* for size */
  3124.  
  3125.       SetDatasetForm(dataset, formptr[n]);
  3126.       RegisterDataset(dataset);
  3127.     } /* outer for */
  3128.  
  3129.   Free(indices);
  3130.   Free(nvar);
  3131.   fclose(inFile);
  3132.   return ngrid;
  3133. } /* RdFmtFun */
  3134.  
  3135. /*********************formatted p3d reading routines, end *********************/
  3136.  
  3137. /*********************unformatted p3d reading routines, start *********************/
  3138. /* only available when there's fortran compiler */
  3139.  
  3140. #ifdef FORTRAN
  3141.  
  3142. #ifdef PROTO
  3143. static int RdUnfmtXYZ(int ndim, int isMGrid, long ***dimsptr,
  3144.                       p3d_grid *gridptr, ObjPtr **formptr, int ***ipp);
  3145. static int RdUnfmtSol(int ndim, int isMGrid, int ngrid, int blank, long **dimsptr,
  3146.               p3d_q *qptr, ObjPtr *formptr, int **ipp);
  3147. static int RdUnfmtFun(int ndim, int isMGrid, int ngrid, int blank, long **dimsptr,
  3148.               p3d_func *funcptr, ObjPtr *formptr, int **ipp);
  3149. #ifdef FORTRAN_
  3150. /* fortran subroutines */
  3151. extern void foropn_(int *lun, int intname[], int *length, int *opstat);
  3152. extern void forcls_(int *lun);
  3153. extern void rngrid_(int *lun, int *ngrid, int *opstat);
  3154. extern void rddims_(int *lun, int *ndim, int *ngrid, int *indims, int *opstat);
  3155. extern void rdgrid_(int *lun, int *ndim, int *indims, int *isiblk, int *iblank,
  3156.             int *iperm, float *gdvals, long *size, int *opstat); 
  3157. extern void rdsolu_(int *lun, int *ndim, int *indims, int *iperm, 
  3158.              float *slvals, long *size, int *opstat);
  3159. extern void rdtime_(int *lun, float *time, int *opstat);
  3160. extern void rdfunc_(int *lun, int *ndim, int *indims, int *nvar, int *iperm,
  3161.             float *fnvals, long *size, int *opstat); 
  3162. #else
  3163. /* fortran subroutines */
  3164. extern void foropn(int *lun, int iname[], int *length, int *opstat);
  3165. extern void forcls(int *lun);
  3166. extern void rngrid(int *lun, int *ngrid, int *opstat);
  3167. extern void rddims(int *lun, int *ndim, int *ngrid, int *indims, int *opstat);
  3168. extern void rdgrid(int *lun, int *ndim, int *indims, int *isiblk, int *iblank,
  3169.            int *iperm, float *gdvals, long *size, int *opstat); 
  3170. extern void rdsolu(int *lun, int *ndim, int *indims, int *iperm, 
  3171.             float *slvals, long *size, int *opstat);
  3172. extern void rdtime(int *lun, float *time, int *opstat);
  3173. extern void rdfunc(int *lun, int *ndim, int *indims, int *nvar, int *iperm,
  3174.            float *fnvals, long *size, int *opstat); 
  3175. #endif
  3176.  
  3177. #else
  3178.  
  3179. static int RdUnfmtXYZ();
  3180. static int RdUnfmtSol();
  3181. static int RdUnfmtFun();
  3182. #ifdef FORTRAN_
  3183. /* fortran subroutines */
  3184. extern void foropn_();
  3185. extern void forcls_();
  3186. extern void rngrid_();
  3187. extern void rddims_();
  3188. extern void rdgrid_();
  3189. extern void rdsolu_();
  3190. extern void rdtime_();
  3191. extern void rdfunc_();
  3192. #else
  3193. /* fortran subroutines */
  3194. extern void foropn();
  3195. extern void forcls();
  3196. extern void rngrid();
  3197. extern void rddims();
  3198. extern void rdgrid();
  3199. extern void rdsolu();
  3200. extern void rdtime();
  3201. extern void rdfunc();
  3202. #endif
  3203. #endif
  3204.  
  3205.  
  3206. /* for reading fortran unformatted grid file */
  3207. #ifdef PROTO
  3208. static int RdUnfmtXYZ(int ndim, int isMGrid, long ***dimsptr,
  3209.                                p3d_grid *gridptr, ObjPtr **formptr, int ***ipp)
  3210. #else
  3211. static int RdUnfmtXYZ(ndim, isMGrid, dimsptr, gridptr, formptr, ipp)
  3212. int ndim;
  3213. int isMGrid;
  3214. long ***dimsptr;
  3215. p3d_grid *gridptr;
  3216. ObjPtr **formptr;
  3217. int ***ipp;
  3218. #endif
  3219. {
  3220.   int ins[MAXNAMELEN]; /* ascii coding of file name */
  3221.   int namelen;         /* file name length */
  3222.   int lun = 2;         /* logical unit number of fortran opened file */
  3223.   int opstat=0;        /* operation status of fortran file op subroutine */
  3224.   int *indims, ngrid;
  3225.   int i, j;
  3226.   
  3227. /* try to open the xyz file, ascii coding of the xyz file name */
  3228.   namelen = (int) strlen(gridptr->name);
  3229.   for (i = 0; i <= namelen; i++)
  3230.      {
  3231.        ins[i] = gridptr->name[i];
  3232.      }
  3233. /* open file with fortran subroutine */
  3234. #ifdef FORTRAN_
  3235.   foropn_(&lun, ins, &namelen, &opstat);
  3236. #else
  3237.   foropn(&lun, ins, &namelen, &opstat);
  3238. #endif
  3239.  
  3240.   if (!opstat)
  3241.     {
  3242.       fprintf(stderr, "Unable to open file %s, aborted\n", gridptr->name);
  3243.       return 0;
  3244.     }
  3245.  
  3246. /* reading ngrid */
  3247.   if (isMGrid)
  3248.     {
  3249.       opstat = 0;
  3250. #ifdef FORTRAN_
  3251.       rngrid_(&lun, &ngrid, &opstat);
  3252. #else
  3253.       rngrid(&lun, &ngrid, &opstat);
  3254. #endif
  3255.       if (opstat <= 0 || ngrid <=0)
  3256.         {
  3257.       fprintf(stderr, "%s format error, invalid NGRID value\n",
  3258.       gridptr->name);
  3259. #ifdef FORTRAN_
  3260.           forcls_(&lun);
  3261. #else
  3262.           forcls(&lun);
  3263. #endif
  3264.       return 0;
  3265.         }
  3266.     }
  3267.   else
  3268.     {
  3269.       ngrid = 1;
  3270.     }
  3271.  
  3272. /* reading dimensions */
  3273.   indims = (int *) 0;
  3274.   *dimsptr = (long **) 0;
  3275.   if (!(indims = (int *) Alloc(ndim*ngrid*sizeof(int))) ||
  3276.       !(*dimsptr = (long **) Alloc(ndim * sizeof(long *))))
  3277.     {
  3278.       fprintf(stderr, "Unable to alloc mem, aborted\n");
  3279. #ifdef FORTRAN_
  3280.           forcls_(&lun);
  3281. #else
  3282.           forcls(&lun);
  3283. #endif
  3284.       Free(indims);
  3285.       return 0;
  3286.     }
  3287.   opstat=0;
  3288. #ifdef FORTRAN_
  3289.   rddims_(&lun, &ndim, &ngrid, indims, &opstat);
  3290. #else
  3291.   rddims(&lun, &ndim, &ngrid, indims, &opstat);
  3292. #endif
  3293.   if (opstat > 0)
  3294.     {
  3295.       for (i = 0; i < ngrid; i++)
  3296.      {
  3297.            if (!((*dimsptr)[i] = (long *) Alloc(ndim * sizeof(long))))
  3298.          {
  3299.            fprintf(stderr, "Unable to alloc memory, aborted\n");
  3300.            i--;
  3301.            opstat = 0;  /* turn off opstat to indicate failure */
  3302.            break;
  3303.              }
  3304.        for (j = 0; j < ndim; j++)
  3305.              {
  3306.            *((*dimsptr)[i]+j)= (long) indims[i*ndim + j];
  3307.              }
  3308.          }
  3309.     }
  3310.   if (opstat <= 0)
  3311.     {
  3312.       fprintf(stderr, "Unable to Read Dimension Values in %s\n",
  3313.           gridptr->name);
  3314. #ifdef FORTRAN_
  3315.           forcls_(&lun);
  3316. #else
  3317.           forcls(&lun);
  3318. #endif
  3319.       freedims(i, *dimsptr);
  3320.       Free(indims);
  3321.       return 0;
  3322.     }
  3323.  
  3324. /* reading grid values*/
  3325.     {
  3326.       long *indices;
  3327.       real *bounds;
  3328.       char *name;   /* dataset name, extracted from file name */
  3329.  
  3330.       /* extracting the dataset name from file name */
  3331.       for(i = namelen -1;
  3332.       i >= 0 && gridptr->name[i] != '/';
  3333.       i--) ;
  3334.       if (i == -1 || gridptr->name[i] == '/')
  3335.     i++;
  3336.       name = gridptr->name + i;
  3337.       for (namelen = 0; 
  3338.        name[namelen] != '\0' && name[namelen] != '.';
  3339.        namelen++) ;
  3340.          
  3341.       /* memory allocation for indices, bounds, formptr, ipp  */
  3342.       indices = (long *) 0;
  3343.       bounds = (real *) 0;
  3344.       *formptr = (ObjPtr *) 0;
  3345.       *ipp = (int **) 0;
  3346.       if ( !(indices = (long *) Alloc(ndim * sizeof(long))) ||
  3347.        !(bounds = (real *) Alloc(ndim * 2 * sizeof(real))) ||
  3348.        !(*formptr = (ObjPtr *) Alloc(ngrid * sizeof(ObjPtr))) ||
  3349.        !(*ipp = (int **) Alloc(ngrid *sizeof(int *))) )
  3350.            /* alloc *ipp regardless of the value of gridptr->blank, but
  3351.           remember to free it if it is of no use */
  3352.         {
  3353.           fprintf(stderr, "Unable to alloc mem, aborted\n");
  3354. #ifdef FORTRAN_
  3355.           forcls_(&lun);
  3356. #else
  3357.           forcls(&lun);
  3358. #endif
  3359.           Free(indims);
  3360.       freedims(ngrid, *dimsptr);
  3361.       Free(*formptr);
  3362.       Free(*ipp);
  3363.       Free(indices);
  3364.       Free(bounds);
  3365.           return 0;
  3366.     }
  3367.  
  3368.       for(i = 0; i < ngrid; i++)
  3369.     {
  3370.           long size, count;
  3371.           float *gdvals;
  3372.       ObjPtr dataset;
  3373.       char datasetname[MAXNAMELEN];
  3374.  
  3375.       /* create a dataset for the grid to be read in */
  3376.           strncpy(datasetname, name, namelen);
  3377.       datasetname[namelen] = '\0';
  3378.       if (ngrid > 1)
  3379.         sprintf(datasetname + namelen, "%d", i + 1);
  3380.           dataset=NewStructuredDataset(datasetname, ndim, (*dimsptr)[i], ndim);
  3381.       if (!SetCurField(FIELD1, dataset))
  3382.         {
  3383.                  fprintf(stderr, "Unable to SETCURFIELD, aborted\n");
  3384. #ifdef FORTRAN_
  3385.                  forcls_(&lun);
  3386. #else
  3387.                  forcls(&lun);
  3388. #endif
  3389.                  Free(indims);
  3390.              freedims(ngrid, *dimsptr);
  3391.          if (gridptr->blank == P3D_BLANK)
  3392.            {
  3393.                  freeiblank(i-1, *ipp);
  3394.            }
  3395.                  else
  3396.            {
  3397.              Free(*ipp);
  3398.            }
  3399.              Free(*formptr);
  3400.              Free(indices);
  3401.              Free(bounds);
  3402.                  return 0;
  3403.         }
  3404.  
  3405.           /* calculate size of one component */
  3406.           size=1;
  3407.           for(j = 0; j < ndim; j++)
  3408.            {
  3409.              size *= (*dimsptr)[i][j];
  3410.            }
  3411.       (*ipp)[i] = (int *) 0;
  3412.           if (gridptr->blank == P3D_BLANK)
  3413.         {
  3414.           if (!((*ipp)[i] = (int *) Alloc(size*sizeof(int))))
  3415.         {
  3416.                  fprintf(stderr, "Unable to alloc mem, aborted\n");
  3417. #ifdef FORTRAN_
  3418.                  forcls_(&lun);
  3419. #else
  3420.                  forcls(&lun);
  3421. #endif
  3422.                  Free(indims);
  3423.              freedims(ngrid, *dimsptr);
  3424.              freeiblank(i, *ipp);
  3425.              Free(indices);
  3426.              Free(bounds);
  3427.              Free(*formptr);
  3428.                  return 0;
  3429.                 }
  3430.         }
  3431.  
  3432.       /* read in grid values */
  3433.           if (!(gdvals = (float *) Alloc(ndim*size*sizeof(float))))
  3434.         {
  3435.                  fprintf(stderr, "Unable to alloc mem, aborted\n");
  3436. #ifdef FORTRAN_
  3437.                  forcls_(&lun);
  3438. #else
  3439.                  forcls(&lun);
  3440. #endif
  3441.                  Free(indims);
  3442.              freedims(ngrid, *dimsptr);
  3443.          if (gridptr->blank == P3D_BLANK)
  3444.            {
  3445.                  freeiblank(i, *ipp);
  3446.            }
  3447.                  else
  3448.            {
  3449.              Free(*ipp);
  3450.            }
  3451.              Free(indices);
  3452.              Free(bounds);
  3453.              Free(*formptr);
  3454.                  return 0;
  3455.         }
  3456.       opstat=0;
  3457. #ifdef FORTRAN_
  3458.           rdgrid_(&lun, &ndim, indims+i*ndim, &gridptr->blank, (*ipp)[i],
  3459.           &(gridptr->permutation), gdvals, &size, &opstat);
  3460. #else
  3461.           rdgrid(&lun, &ndim, indims+i*ndim, &gridptr->blank, (*ipp)[i],
  3462.           &(gridptr->permutation), gdvals, &size, &opstat);
  3463. #endif
  3464.           if (!opstat)
  3465.         {
  3466.           fprintf(stderr, "Read error on %d grid of file %s\n", 
  3467.                i, gridptr->name);
  3468. #ifdef FORTRAN_
  3469.               forcls_(&lun);
  3470. #else
  3471.               forcls(&lun);
  3472. #endif
  3473.               Free(indims);
  3474.           freedims(ngrid, *dimsptr);
  3475.           if (gridptr->blank == P3D_BLANK)
  3476.             {
  3477.               freeiblank(i, *ipp);
  3478.         }
  3479.               else
  3480.             {
  3481.           Free(*ipp);
  3482.         }
  3483.           Free(indices);
  3484.           Free(bounds);
  3485.           Free(*formptr);
  3486.           Free(gdvals);
  3487.               return 0;
  3488.         }
  3489.  
  3490.  
  3491.       /* assign gdvals to dataset and calculate bounds here */
  3492.           for(j = 0; j < ndim; j++)
  3493.         {
  3494.           indices[j] = 0;
  3495.           /* give bounds initial values */
  3496.           bounds[j*2] = bounds[j*2 + 1] = gdvals[j*size];
  3497.             }
  3498.           for(count = 0; count < size; count++)
  3499.         {
  3500.           int k, whichcomponent;
  3501.  
  3502.           /* assign grid values */
  3503.           for(whichcomponent = 0; whichcomponent < ndim; whichcomponent++)
  3504.         {
  3505.           real rval;
  3506.  
  3507.           rval = (real) gdvals[whichcomponent*size + count];
  3508.           PutFieldComponent(FIELD1, whichcomponent, indices, rval);
  3509.                   if (rval < bounds[whichcomponent*2])
  3510.             bounds[whichcomponent*2] = rval;
  3511.                   else if (rval > bounds[whichcomponent*2 + 1])
  3512.             bounds[whichcomponent*2 + 1] = rval;
  3513.         }
  3514.           /* update indices */
  3515.           /* row major */
  3516.               for(k = ndim - 1; k >= 0; k--)
  3517.         {
  3518.           indices[k]++;
  3519.           if (indices[k] == (*dimsptr)[i][k])
  3520.             {
  3521.               indices[k] = 0;
  3522.               continue;
  3523.             }
  3524.           break;
  3525.         }
  3526.         }
  3527.  
  3528.           Free(gdvals);
  3529.       (*formptr)[i]=NewCurvilinearDataForm(datasetname, ndim,
  3530.                                (*dimsptr)[i], bounds, dataset);
  3531.     }
  3532.       Free(indices);
  3533.       Free(bounds);
  3534.     }
  3535.  
  3536. /* close file with fortran subroutine */
  3537. #ifdef FORTRAN_
  3538.   forcls_(&lun);
  3539. #else
  3540.   forcls(&lun);
  3541. #endif
  3542.  
  3543.   Free(indims); 
  3544.   if (gridptr->blank != P3D_BLANK)
  3545.     Free(*ipp);
  3546.   return i;
  3547. } /* RdUnfmtXYZ */
  3548.  
  3549.  
  3550. /* for reading fortran unformatted solution file */
  3551. #ifdef PROTO
  3552. static int RdUnfmtSol(int ndim, int isMGrid, int ngrid, int blank, long **dimsptr,
  3553.               p3d_q *qptr, ObjPtr *formptr, int **ipp)
  3554. #else
  3555. static int RdUnfmtSol(ndim, isMGrid, ngrid, blank, dimsptr, qptr, formptr, ipp)
  3556. int ndim;
  3557. int isMGrid;
  3558. int ngrid;
  3559. int blank;
  3560. long **dimsptr;
  3561. p3d_q *qptr;
  3562. ObjPtr *formptr;
  3563. int **ipp;
  3564. #endif
  3565. {
  3566.   int ins[MAXNAMELEN]; /* ascii coding of file name */
  3567.   int namelen;         /* file name length */
  3568.   int lun = 2;         /* logical unit number of fortran opened file */
  3569.   int opstat=0;        /* operation status of fortran file op subroutine */
  3570.   int *indims;
  3571.   int numsets;         /* number of datasets read so far */
  3572.   int i, j;
  3573.   int myngrid;         /* the ngrid read from the solution file */
  3574.   long *indices;
  3575.   char *name;   /* dataset name, extracted from file name */
  3576.  
  3577.  
  3578. /* try to open the solution file, ascii coding of the solution file name */
  3579.   namelen = (int) strlen(qptr->name);
  3580.   for (i = 0; i <= namelen; i++)
  3581.      {
  3582.        ins[i] = qptr->name[i];
  3583.      }
  3584. /* open file with fortran subroutine */
  3585. #ifdef FORTRAN_
  3586.   foropn_(&lun, ins, &namelen, &opstat);
  3587. #else
  3588.   foropn(&lun, ins, &namelen, &opstat);
  3589. #endif
  3590.   if (!opstat)
  3591.     {
  3592.       fprintf(stderr, "Unable to open file %s, aborted\n", qptr->name);
  3593.       return 0;
  3594.     }
  3595.  
  3596.   /* extracting the dataset name from file name */
  3597.   for(i = namelen -1; i >= 0 && qptr->name[i] != '/'; i--)
  3598.     ;
  3599.   if (i == -1 || qptr->name[i] == '/')
  3600.     i++;
  3601.   name = qptr->name + i;
  3602.  
  3603.   for (namelen = 0; 
  3604.        name[namelen] != '\0' && name[namelen] != '.';
  3605.        namelen++) ;
  3606.  
  3607.   /* allocate storage for indices */
  3608.   if (!(indices = (long *) Alloc(ndim*sizeof(long))))
  3609.     {
  3610.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  3611. #ifdef FORTRAN_
  3612.       forcls_(&lun);
  3613. #else
  3614.       forcls(&lun);
  3615. #endif
  3616.       return 0;
  3617.     }
  3618.  
  3619.   /* 0 dataset read so far */
  3620.   numsets = 0;
  3621.  
  3622. /* do multiple datasets starting here */
  3623.   do
  3624.     {
  3625.       /* reading myngrid */
  3626.       if (isMGrid)
  3627.         {
  3628.           opstat = 0;
  3629. #ifdef FORTRAN_
  3630.           rngrid_(&lun, &myngrid, &opstat);
  3631. #else
  3632.           rngrid(&lun, &myngrid, &opstat);
  3633. #endif
  3634.           if (opstat == EOF)
  3635.         {
  3636.           if (numsets == 0)
  3637.             {
  3638.               fprintf(stderr, "Error in %s: Empty file\n", qptr->name);
  3639.             }
  3640. #ifdef FORTRAN_
  3641.               forcls_(&lun);
  3642. #else
  3643.               forcls(&lun);
  3644. #endif
  3645.               Free(indices);
  3646.           /* if numsets >= 1 then this is an EOF, else error */
  3647.           return (numsets ? ngrid : 0);
  3648.         }
  3649.           else if (!opstat || myngrid <=0)
  3650.             {
  3651.           fprintf(stderr, "%s format error, invalid NGRID value\n", qptr->name);
  3652. #ifdef FORTRAN_
  3653.               forcls_(&lun);
  3654. #else
  3655.               forcls(&lun);
  3656. #endif
  3657.               Free(indices);
  3658.           return 0;
  3659.             }
  3660.         }
  3661.       else
  3662.         {
  3663.           myngrid = 1;
  3664.         }
  3665.     
  3666.       /* if myngrid is not the same as ngrid then error */
  3667.       if (myngrid != ngrid)
  3668.         {
  3669.           fprintf(stderr, "Error in %s: NGRID value does not match that"
  3670.               "of its XYZ file\n", qptr->name);
  3671. #ifdef FORTRAN_
  3672.           forcls_(&lun);
  3673. #else
  3674.           forcls(&lun);
  3675. #endif
  3676.           Free(indices);
  3677.           return 0;
  3678.         }
  3679.     
  3680.        /* reading dimensions */
  3681.       indims = (int *) 0;
  3682.       if (!(indims = (int *) Alloc(ndim*myngrid*sizeof(int))))
  3683.         {
  3684.           fprintf(stderr, "Unable to alloc mem, aborted\n");
  3685. #ifdef FORTRAN_
  3686.           forcls_(&lun);
  3687. #else
  3688.           forcls(&lun);
  3689. #endif
  3690.           Free(indices);
  3691.           return 0;
  3692.         }
  3693.       opstat=0;
  3694. #ifdef FORTRAN_
  3695.       rddims_(&lun, &ndim, &ngrid, indims, &opstat);
  3696. #else
  3697.       rddims(&lun, &ndim, &ngrid, indims, &opstat);
  3698. #endif
  3699.       if (opstat == EOF && !isMGrid)
  3700.         {
  3701.           if (numsets == 0)
  3702.         {
  3703.           fprintf(stderr, "Error in %s: Empty file\n", qptr->name);
  3704.         }
  3705. #ifdef FORTRAN_
  3706.           forcls_(&lun);
  3707. #else
  3708.           forcls(&lun);
  3709. #endif
  3710.           Free(indices);
  3711.           Free(indims);
  3712.           return (numsets ? ngrid : 0);
  3713.         }
  3714.       else if (opstat > 0)
  3715.         {
  3716.           /* Check dimension consistency with dims read from XYZ file */
  3717.           for (i = 0; i < ngrid; i++)
  3718.          {
  3719.            for (j = 0; j < ndim; j++)
  3720.                  {
  3721.                if (dimsptr[i][j] != (long) indims[i*ndim + j])
  3722.              {
  3723.                fprintf(stderr, "Dimensions mismatch between %s "
  3724.                "and its XYZ file\n", qptr->name);
  3725. #ifdef FORTRAN_
  3726.                        forcls_(&lun);
  3727. #else
  3728.                        forcls(&lun);
  3729. #endif
  3730.                        Free(indices);
  3731.                        Free(indims);
  3732.                        return 0;
  3733.                      }
  3734.                  }
  3735.              }
  3736.         }
  3737.       else
  3738.       /* opstat is 0 or (opstat is EOF and isMGrid is true) */
  3739.         {
  3740.           fprintf(stderr, "Unable to Read Dimension Values in %s\n", qptr->name);
  3741. #ifdef FORTRAN_
  3742.               forcls_(&lun);
  3743. #else
  3744.               forcls(&lun);
  3745. #endif
  3746.           Free(indices);
  3747.           Free(indims);
  3748.           return 0;
  3749.         }
  3750.  
  3751.       for(i = 0; i < ngrid; i++)
  3752.     {
  3753.           long size, count;
  3754.           float *slvals;
  3755.       ObjPtr density, pressure, dataset;
  3756.       char densityname[MAXNAMELEN];
  3757.       char pressurename[MAXNAMELEN];
  3758.       char datasetname[MAXNAMELEN];
  3759.       float time;
  3760.       int whichcomponent;
  3761.  
  3762.       /* get the time data of the current grid solutions */
  3763.           opstat = 0;
  3764. #ifdef FORTRAN_
  3765.           rdtime_(&lun, &time, &opstat);
  3766. #else
  3767.           rdtime(&lun, &time, &opstat);
  3768. #endif
  3769.           if (!opstat)
  3770.         {
  3771.           fprintf(stderr, "Error in %s, invalid FAMACH...TIME value\n",
  3772.               qptr->name);
  3773. #ifdef FORTRAN_
  3774.               forcls_(&lun);
  3775. #else
  3776.               forcls(&lun);
  3777. #endif
  3778.               Free(indices);
  3779.               Free(indims);
  3780.               return 0;
  3781.         }
  3782.  
  3783.           /* calculate size of one component */
  3784.           size=1;
  3785.           for(j = 0; j < ndim; j++)
  3786.            {
  3787.              size *= dimsptr[i][j];
  3788.            }
  3789.       /* read in solution, density, and pressure values */
  3790.           if (!(slvals = (float *) Alloc((ndim+2)*size*sizeof(float))))
  3791.         {
  3792.                  fprintf(stderr, "Unable to alloc mem, aborted\n");
  3793. #ifdef FORTRAN_
  3794.                  forcls_(&lun);
  3795. #else
  3796.                  forcls(&lun);
  3797. #endif
  3798.                  Free(indims);
  3799.              Free(indices);
  3800.                  return 0;
  3801.         }
  3802.       opstat=0;
  3803. #ifdef FORTRAN_
  3804.           rdsolu_(&lun, &ndim, indims+i*ndim, &(qptr->permutation), 
  3805.           slvals, &size, &opstat);
  3806. #else
  3807.           rdsolu(&lun, &ndim, indims+i*ndim, &(qptr->permutation), 
  3808.           slvals, &size, &opstat);
  3809. #endif
  3810.           if (!opstat)
  3811.         {
  3812.           fprintf(stderr, "Read error on file %s\n", qptr->name);
  3813. #ifdef FORTRAN_
  3814.               forcls_(&lun);
  3815. #else
  3816.               forcls(&lun);
  3817. #endif
  3818.               Free(indims);
  3819.           Free(indices);
  3820.           Free(slvals);
  3821.               return 0;
  3822.         }
  3823.  
  3824.       /* create datasets for solution, density, pressure */
  3825.           strncpy(datasetname, name, namelen);
  3826.       datasetname[namelen] = '\0';
  3827.       if (ngrid > 1)
  3828.         sprintf(datasetname + namelen, "%d", i + 1);
  3829.           sprintf(densityname, "%s%s", datasetname, " Density");
  3830.           sprintf(pressurename, "%s%s", datasetname, " Pressure");
  3831.           density=NewStructuredDataset(densityname, ndim, dimsptr[i], 0);
  3832.           pressure=NewStructuredDataset(pressurename, ndim, dimsptr[i], 0);
  3833.           dataset=NewStructuredDataset(datasetname, ndim, dimsptr[i], ndim);
  3834.  
  3835.       /* assign slvals to density, pressure, and dataset */
  3836.           for(whichcomponent = 0; whichcomponent < ndim + 2; whichcomponent++)
  3837.         {
  3838.           ObjPtr whichdataset;
  3839.  
  3840.           if (whichcomponent == 0)
  3841.         /* for density */
  3842.         {
  3843.           whichdataset = density;
  3844.             }
  3845.               else if (whichcomponent == 1)
  3846.         /* for pressure */
  3847.         {
  3848.           whichdataset = pressure;
  3849.         }
  3850.               else
  3851.             /* for solution */
  3852.         {
  3853.           whichdataset = dataset;
  3854.         }
  3855.  
  3856.           if (!SetCurField(FIELD1, whichdataset))
  3857.             {
  3858.                   fprintf(stderr, "Unable to SETCURFIELD, aborted\n");
  3859. #ifdef FORTRAN_
  3860.                   forcls_(&lun);
  3861. #else
  3862.                   forcls(&lun);
  3863. #endif
  3864.                   Free(indims);
  3865.               Free(indices);
  3866.           Free(slvals);
  3867.                   return 0;
  3868.             }
  3869.               for(j = 0; j < ndim; j++)
  3870.             {
  3871.               indices[j] = 0;
  3872.                 }
  3873.               for(count = 0; count < size; count++)
  3874.             {
  3875.               int k, iblank = 1;
  3876.           real rval;
  3877.  
  3878.               rval = (real) slvals[whichcomponent*size + count];
  3879.               /* do iblank if needed */
  3880.               if (blank == P3D_BLANK)
  3881.                 {
  3882.                   getIblank(&iblank, ipp, indices, dimsptr[i], ndim, i);
  3883.             }
  3884.               /* do check of negative values of density or pressure 
  3885.              if required */
  3886.           if (whichcomponent < 2 && qptr->check == P3D_CHECK && 
  3887.                rval < (real) 0.0)
  3888.             {
  3889.               fprintf(stderr, "Negative %s value %f found in Q file %s and is converted to %f\n", whichdataset == pressure ? "pressure" : "density", rval, qptr->name, 0.0);
  3890.               rval = (real) 0.0;
  3891.             }
  3892.           PutFieldComponent(FIELD1, 
  3893.              whichcomponent < 2 ? 0 : whichcomponent - 2, indices, 
  3894.              iblank == 0 ? (real) missingData : rval);
  3895.               /* update indices */
  3896.           /* row major way */
  3897.                   for(k = ndim - 1; k >= 0; k--)
  3898.             {
  3899.               indices[k]++;
  3900.               if (indices[k] == dimsptr[i][k])
  3901.                 {
  3902.                   indices[k] = 0;
  3903.                   continue;
  3904.                 }
  3905.               break;
  3906.             }
  3907.         } /* for count = 0 ~ size */
  3908.         } /* for whichcomponent = 0 ~ ndim+2 */
  3909.  
  3910.           Free(slvals);
  3911.       SetDatasetForm(density, formptr[i]);
  3912.       SetDatasetForm(pressure, formptr[i]);
  3913.       SetDatasetForm(dataset, formptr[i]);
  3914.       if ( qptr->isMDataset )
  3915.         {
  3916.           RegisterTimeDataset(density, (real) time);
  3917.           RegisterTimeDataset(pressure, (real) time);
  3918.           RegisterTimeDataset(dataset, (real) time);
  3919.         }
  3920.           else
  3921.         {
  3922.           RegisterDataset(density);
  3923.           RegisterDataset(pressure);
  3924.           RegisterDataset(dataset);
  3925.         }
  3926.     }/* for loop for i = 0, ngrid - 1 */
  3927.       numsets++;
  3928.     } while(qptr->isMDataset);
  3929.  
  3930. /* close file with fortran subroutine */
  3931. #ifdef FORTRAN_
  3932.   forcls_(&lun);
  3933. #else
  3934.   forcls(&lun);
  3935. #endif
  3936.  
  3937.   Free(indims); 
  3938.   Free(indices);
  3939.   return ngrid;
  3940. } /* RdUnfmtSol */
  3941.  
  3942. /* for reading fortran unformatted function file */
  3943. #ifdef PROTO
  3944. static int RdUnfmtFun(int ndim, int isMGrid, int ngrid, int blank, long **dimsptr,
  3945.               p3d_func *funcptr, ObjPtr *formptr, int **ipp)
  3946. #else
  3947. static int RdUnfmtFun(ndim, isMGrid, ngrid, blank, dimsptr, funcptr, formptr, ipp)
  3948. int ndim;
  3949. int isMGrid;
  3950. int ngrid;
  3951. int blank;
  3952. long **dimsptr;
  3953. p3d_func *funcptr;
  3954. ObjPtr *formptr;
  3955. int **ipp;
  3956. #endif
  3957. {
  3958.   int ins[MAXNAMELEN]; /* ascii coding of file name */
  3959.   int namelen;         /* file name length */
  3960.   int lun = 2;         /* logical unit number of fortran opened file */
  3961.   int opstat=0;        /* operation status of fortran file op subroutine */
  3962.   int *indims;
  3963.   int i, j;
  3964.   int myngrid;         /* the ngrid read from the solution file */
  3965.   int nvardim;         /* number of dimensions + 1 */
  3966.   long *indices;
  3967.   char *name;   /* dataset name, extracted from file name */
  3968.  
  3969.  
  3970. /* try to open the solution file, ascii coding of the solution file name */
  3971.   namelen = (int) strlen(funcptr->name);
  3972.   for (i = 0; i <= namelen; i++)
  3973.      {
  3974.        ins[i] = funcptr->name[i];
  3975.      }
  3976. /* open file with fortran subroutine */
  3977. #ifdef FORTRAN_
  3978.   foropn_(&lun, ins, &namelen, &opstat);
  3979. #else
  3980.   foropn(&lun, ins, &namelen, &opstat);
  3981. #endif
  3982.   if (!opstat)
  3983.     {
  3984.       fprintf(stderr, "Unable to open file %s, aborted\n", funcptr->name);
  3985.       return 0;
  3986.     }
  3987.  
  3988.   /* read myngrid */
  3989.   if (isMGrid)
  3990.     {
  3991.       opstat = 0;
  3992. #ifdef FORTRAN_
  3993.       rngrid_(&lun, &myngrid, &opstat);
  3994. #else
  3995.       rngrid(&lun, &myngrid, &opstat);
  3996. #endif
  3997.       if (opstat <= 0 || myngrid <=0)
  3998.         {
  3999.       fprintf(stderr, "%s format error, invalid NGRID value\n", funcptr->name);
  4000. #ifdef FORTRAN_
  4001.           forcls_(&lun);
  4002. #else
  4003.           forcls(&lun);
  4004. #endif
  4005.       return 0;
  4006.         }
  4007.     }
  4008.   else
  4009.     {
  4010.       myngrid = 1;
  4011.     }
  4012.  
  4013.   /* if myngrid is not the same as ngrid then error */
  4014.   if (myngrid != ngrid)
  4015.     {
  4016.       fprintf(stderr, "Error in %s: NGRID value does not match that"
  4017.           "of its XYZ file\n", funcptr->name);
  4018. #ifdef FORTRAN_
  4019.       forcls_(&lun);
  4020. #else
  4021.       forcls(&lun);
  4022. #endif
  4023.       return 0;
  4024.     }
  4025.  
  4026.   /* extracting the dataset name from file name */
  4027.   for(i = namelen -1; i >= 0 && funcptr->name[i] != '/'; i--)
  4028.     ;
  4029.   if (i == -1 || funcptr->name[i] == '/')
  4030.     i++;
  4031.   name = funcptr->name + i;
  4032.  
  4033.   for (namelen = 0; 
  4034.        name[namelen] != '\0' && name[namelen] != '.';
  4035.        namelen++) ;
  4036.  
  4037.   /* allocate storage for indices */
  4038.   if (!(indices = (long *) Alloc(ndim*sizeof(long))))
  4039.     {
  4040.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  4041. #ifdef FORTRAN_
  4042.       forcls_(&lun);
  4043. #else
  4044.       forcls(&lun);
  4045. #endif
  4046.       return 0;
  4047.     }
  4048.  
  4049.   /* reading dimensions */
  4050.   indims = (int *) 0;
  4051.   nvardim = ndim + 1;
  4052.   if (!(indims = (int *) Alloc(nvardim*myngrid*sizeof(int))))
  4053.     {
  4054.       fprintf(stderr, "Unable to alloc mem, aborted\n");
  4055. #ifdef FORTRAN_
  4056.       forcls_(&lun);
  4057. #else
  4058.       forcls(&lun);
  4059. #endif
  4060.       Free(indices);
  4061.       return 0;
  4062.     }
  4063.   opstat=0;
  4064. #ifdef FORTRAN_
  4065.   rddims_(&lun, &nvardim, &ngrid, indims, &opstat);
  4066. #else
  4067.   rddims(&lun, &nvardim, &ngrid, indims, &opstat);
  4068. #endif
  4069.   if (opstat > 0)
  4070.     {
  4071.       /* Check dimension consistency with dims read from XYZ file */
  4072.       for (i = 0; i < ngrid; i++)
  4073.      {
  4074.        for (j = 0; j < ndim; j++)
  4075.              {
  4076.            if (dimsptr[i][j] != (long) indims[i*nvardim + j])
  4077.          {
  4078.            fprintf(stderr, "Dimensions mismatch between %s "
  4079.            "and its XYZ file\n", funcptr->name);
  4080. #ifdef FORTRAN_
  4081.                    forcls_(&lun);
  4082. #else
  4083.                    forcls(&lun);
  4084. #endif
  4085.                    Free(indices);
  4086.                    Free(indims);
  4087.                    return 0;
  4088.                  }
  4089.              }
  4090.          }
  4091.     }
  4092.   else
  4093.   /* opstat <= 0 */
  4094.     {
  4095.       fprintf(stderr, "Unable to Read Dimension Values in %s\n", funcptr->name);
  4096. #ifdef FORTRAN_
  4097.           forcls_(&lun);
  4098. #else
  4099.           forcls(&lun);
  4100. #endif
  4101.       Free(indices);
  4102.       Free(indims);
  4103.       return 0;
  4104.     }
  4105.  
  4106.       for(i = 0; i < ngrid; i++)
  4107.     {
  4108.           long size, count;
  4109.           float *fnvals;
  4110.       ObjPtr dataset;
  4111.       char datasetname[MAXNAMELEN];
  4112.  
  4113.           /* calculate size of one component */
  4114.           size=1;
  4115.           for(j = 0; j < ndim; j++)
  4116.            {
  4117.              size *= dimsptr[i][j];
  4118.            }
  4119.       /* read in function values */
  4120.           if (!(fnvals = (float *) Alloc(nvardim*size*sizeof(float))))
  4121.         {
  4122.                  fprintf(stderr, "Unable to alloc mem, aborted\n");
  4123. #ifdef FORTRAN_
  4124.                  forcls_(&lun);
  4125. #else
  4126.                  forcls(&lun);
  4127. #endif
  4128.                  Free(indims);
  4129.              Free(indices);
  4130.                  return 0;
  4131.         }
  4132.       opstat=0;
  4133. #ifdef FORTRAN_
  4134.           rdfunc_(&lun, &ndim, indims+i*nvardim, indims+i*nvardim+ndim, 
  4135.                &(funcptr->permutation), fnvals, &size, &opstat);
  4136. #else
  4137.           rdfunc(&lun, &ndim, indims+i*nvardim, indims+i*nvardim+ndim, 
  4138.                &(funcptr->permutation), fnvals, &size, &opstat);
  4139. #endif
  4140.           if (!opstat)
  4141.         {
  4142.           fprintf(stderr, "Read error on file %s\n", funcptr->name);
  4143. #ifdef FORTRAN_
  4144.               forcls_(&lun);
  4145. #else
  4146.               forcls(&lun);
  4147. #endif
  4148.               Free(indims);
  4149.           Free(indices);
  4150.           Free(fnvals);
  4151.               return 0;
  4152.         }
  4153.  
  4154.       /* create datasets for function */
  4155.           strncpy(datasetname, name, namelen);
  4156.       datasetname[namelen] = '\0';
  4157.       if (ngrid > 1)
  4158.         sprintf(datasetname + namelen, "%d", i + 1);
  4159.           dataset=NewStructuredDataset(datasetname, ndim, dimsptr[i], 
  4160.            indims[i*nvardim+ndim] == 1 ? 0 : indims[i*nvardim+ndim]);
  4161.  
  4162.       /* assign fnvals to dataset */
  4163.           for(j = 0; j < ndim; j++)
  4164.         {
  4165.           indices[j] = 0;
  4166.             }
  4167.           if (!SetCurField(FIELD1, dataset))
  4168.             {
  4169.               fprintf(stderr, "Unable to SETCURFIELD, aborted\n");
  4170. #ifdef FORTRAN_
  4171.               forcls_(&lun);
  4172. #else
  4173.               forcls(&lun);
  4174. #endif
  4175.               Free(indims);
  4176.               Free(indices);
  4177.               Free(fnvals);
  4178.               return 0;
  4179.         }
  4180.           for(count = 0; count < size; count++)
  4181.         {
  4182.           int whichcomponent, k;
  4183.           int iblank = 1;
  4184.  
  4185.           /* do iblank if needed */
  4186.           if (blank == P3D_BLANK)
  4187.             {
  4188.               getIblank(&iblank, ipp, indices, dimsptr[i], ndim, i);
  4189.         }
  4190.           /* assign solution values */
  4191.           for(whichcomponent = 0; whichcomponent < indims[i*nvardim+ndim];
  4192.           whichcomponent++)
  4193.         {
  4194.           real rval;
  4195.  
  4196.           rval = (real) fnvals[whichcomponent*size + count];
  4197.  
  4198.           PutFieldComponent(FIELD1, whichcomponent, indices, 
  4199.           iblank == 0 ? (real) missingData : rval);
  4200.         }
  4201.           /* update indices */
  4202.           /* row major way */
  4203.               for(k = ndim - 1; k >= 0; k--)
  4204.         {
  4205.           indices[k]++;
  4206.           if (indices[k] == dimsptr[i][k])
  4207.             {
  4208.               indices[k] = 0;
  4209.               continue;
  4210.             }
  4211.           break;
  4212.         }
  4213.             }
  4214.           Free(fnvals);
  4215.           SetDatasetForm(dataset, formptr[i]);
  4216.              RegisterDataset(dataset);
  4217.     }/* for loop for i = 0, ngrid - 1 */
  4218.  
  4219. /* close file with fortran subroutine */
  4220. #ifdef FORTRAN_
  4221.   forcls_(&lun);
  4222. #else
  4223.   forcls(&lun);
  4224. #endif
  4225.  
  4226.   Free(indims); 
  4227.   Free(indices);
  4228.   return ngrid;
  4229. } /* RdUnfmtFun */
  4230.  
  4231. #endif /* for ifdef FORTRAN */
  4232.  
  4233. /*********************unformatted p3d reading routines, end ***********************/
  4234.  
  4235. /*********************binary p3d reading routines, start *********************/
  4236.  
  4237. #ifdef PROTO
  4238. static int RdBinXYZ(int ndim, int isMGrid, long ***dimsptr,
  4239.                                p3d_grid *gridptr, ObjPtr **formptr, int ***ipp)
  4240. #else
  4241. static int RdBinXYZ(ndim, isMGrid, dimsptr, gridptr, formptr, ipp)
  4242. int ndim;
  4243. int isMGrid;
  4244. long ***dimsptr;
  4245. p3d_grid *gridptr;
  4246. ObjPtr **formptr;
  4247. int ***ipp;
  4248. #endif
  4249. {
  4250.   FILE *inFile;
  4251.   int ngrid;            /* number of grids */
  4252.  
  4253.   int i, j;
  4254.   int n;
  4255.  
  4256.   char *name;           /* dataset name starting point */
  4257.   int namelen;
  4258.  
  4259.   int ival;             /* an integer value */
  4260.   float fval;           /* a float value */
  4261.  
  4262.   long *indices;
  4263.   real *bounds;
  4264.  
  4265.   if (!(inFile = fopen(gridptr->name, "r")))
  4266.     {
  4267.       fprintf(stderr, "Unable to open %s, aborted\n", gridptr->name);
  4268.       return 0;
  4269.     }
  4270.  
  4271.   /* determine the dataset name now */
  4272.   for (i = (int) strlen(gridptr->name) - 1; i >= 0 && gridptr->name[i] != '/'; i--)
  4273.      ;
  4274.   if (i == -1 || gridptr->name[i] == '/')
  4275.     i++;
  4276.   name = gridptr->name + i;
  4277.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  4278.      ;
  4279.  
  4280.   if (isMGrid)
  4281.     {
  4282.       if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4283.         {
  4284.           fprintf(stderr, "Error reading NGRID value in %s\n", gridptr->name);
  4285.           fclose(inFile);
  4286.           return 0;
  4287.         } 
  4288.       else if (ival <= 0)
  4289.         {
  4290.           fprintf(stderr, "%s format error, invalid NGRID value\n", gridptr->name);
  4291.           fclose(inFile);
  4292.           return 0;
  4293.         } 
  4294.       else
  4295.     {
  4296.       ngrid = ival;
  4297.     }
  4298.     }
  4299.   else
  4300.     {
  4301.       ngrid = 1;
  4302.     }
  4303.   
  4304.   /* alloc dimsptr and indices storage */
  4305.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  4306.     {
  4307.       fprintf(stderr, "Unable to alloc memory\n");
  4308.       fclose(inFile);
  4309.       return 0;
  4310.     }
  4311.   if (!(*dimsptr = (long **) Alloc(ngrid * sizeof(long *))))
  4312.     {
  4313.       fprintf(stderr, "Unable to alloc memory\n");
  4314.       Free(indices);
  4315.       fclose(inFile);
  4316.       return 0;
  4317.     }
  4318.  
  4319.   /* read in all the dimensions */
  4320.   for (i = 0; i < ngrid; i++)
  4321.     {
  4322.       if (!((*dimsptr)[i] = (long *) Alloc(ndim * sizeof(long))))
  4323.         {
  4324.           fprintf(stderr, "Unable to alloc memory\n");
  4325.           Free(indices);
  4326.           freedims(i, *dimsptr);
  4327.           fclose(inFile);
  4328.           return 0;
  4329.         }
  4330.  
  4331.       for (j = 0; j < ndim; j++)
  4332.          {
  4333.            if (1 != fread((void *) &ival, sizeof(int), 1, inFile) || ival <= 0)
  4334.              {
  4335.                fprintf(stderr, "Invalid dimension value in %s\n", gridptr->name);
  4336.                Free(indices);
  4337.                freedims(i, *dimsptr);
  4338.                fclose(inFile);
  4339.                return 0;
  4340.              }
  4341.            *((*dimsptr)[i] + j) = (long) ival; 
  4342.          }
  4343.     }
  4344.  
  4345.   /* alloc strorage for pointer to bounds */
  4346.   if (!(bounds = (real *) Alloc(ndim * 2 * sizeof(real))))
  4347.     {
  4348.       fprintf(stderr, "Unable to alloc memory\n");
  4349.       Free(indices);
  4350.       freedims(ngrid, *dimsptr);
  4351.       fclose(inFile);
  4352.       return 0;
  4353.     }
  4354.   /* alloc storage for pointers to dataforms */
  4355.   if (!(*formptr = (ObjPtr *) Alloc(ngrid * sizeof(ObjPtr))))
  4356.     {
  4357.       fprintf(stderr, "Unable to alloc memory\n");
  4358.       Free(indices);
  4359.       Free(bounds);
  4360.       freedims(ngrid, *dimsptr);
  4361.       fclose(inFile);
  4362.       return 0;
  4363.     }
  4364.  
  4365.   /* alloc iblankptr if IBLANK */
  4366.   if (gridptr->blank == P3D_BLANK)
  4367.     {
  4368.       if (!(*ipp = (int **) Alloc(ngrid * sizeof(int *))))
  4369.         {
  4370.           fprintf(stderr, "Unable to alloc memory\n");
  4371.           Free(indices);
  4372.           Free(bounds);
  4373.           Free(*formptr);
  4374.           freedims(ngrid, *dimsptr);
  4375.           fclose(inFile);
  4376.           return 0;
  4377.         }
  4378.     }
  4379.  
  4380.   for (n = 0; n < ngrid; n++)
  4381.     {
  4382.       long size, count; /* number of grid values promised */
  4383.       ObjPtr dataset;
  4384.       char datasetname[MAXNAMELEN];
  4385.       int whichcomponent;
  4386.  
  4387.       strncpy(datasetname, name, namelen);
  4388.       datasetname[namelen] = '\0';
  4389.       if (ngrid > 1)
  4390.         {
  4391.           sprintf(datasetname + namelen, "%d", n + 1);
  4392.         }
  4393.  
  4394.       /* setup a dataset for the grid data */
  4395.       dataset = NewStructuredDataset(datasetname, ndim, (*dimsptr)[n], ndim);
  4396.       if (!SetCurField(FIELD1, dataset))
  4397.         {
  4398.           fprintf(stderr, "Unable to SetCurField\n");
  4399.           Free(indices);
  4400.           Free(bounds);
  4401.           Free(*formptr);
  4402.           freedims(ngrid, *dimsptr);
  4403.           if (gridptr->blank == P3D_BLANK)
  4404.             {
  4405.               freeiblank(n-1, *ipp);
  4406.             }
  4407.           fclose(inFile);
  4408.           return 0;
  4409.         }
  4410.  
  4411.       size = 1;
  4412.       for (i = 0; i < ndim; i++)
  4413.          {
  4414.            size *= (*dimsptr)[n][i];
  4415.          }
  4416.       size *= (long)(ndim + (gridptr->blank == P3D_BLANK));
  4417.  
  4418.       if (gridptr->blank == P3D_BLANK)
  4419.         {
  4420.           /* alloc storage for the iblank arrays */
  4421.           if (!((*ipp)[n] = (int *) Alloc(size * sizeof(int))))
  4422.             {
  4423.               fprintf(stderr, "Unable to alloc memory\n");
  4424.               freeiblank(n, *ipp);
  4425.               Free(indices);
  4426.               Free(bounds);
  4427.               Free(*formptr);
  4428.               freedims(ngrid, *dimsptr);
  4429.               fclose(inFile);
  4430.               return 0;
  4431.             }
  4432.         }
  4433.  
  4434.       for (i = 0; i < ndim; i++)
  4435.          {
  4436.            indices[i] = 0;
  4437.          }
  4438.       whichcomponent = 0;
  4439.       /* read in all the grid values */
  4440.       for (count = 0; count < size; count++)
  4441.         {
  4442.           int first = 1;
  4443.  
  4444.           if (whichcomponent < ndim) /* read grid values */
  4445.             {
  4446.               if (1 != fread((void *) &fval, sizeof(float), 1, inFile))
  4447.                 {
  4448.               if (feof(inFile))
  4449.             {
  4450.                       fprintf(stderr, "Short of data in %s\n", gridptr->name);
  4451.             }
  4452.                   else
  4453.             {
  4454.                       fprintf(stderr, "Error reading %s\n", gridptr->name);
  4455.             }
  4456.                   Free(indices);
  4457.                   Free(bounds);
  4458.                   Free(*formptr);
  4459.                   freedims(ngrid, *dimsptr);
  4460.                   if (gridptr->blank == P3D_BLANK)
  4461.                     {
  4462.                       freeiblank(n + 1, *ipp);
  4463.                     }
  4464.                   fclose(inFile);
  4465.                   return 0;
  4466.                 }
  4467.               PutFieldComponent(FIELD1, whichcomponent, indices, (real) fval);  
  4468.               for (j = 0; first && j < ndim; j++)
  4469.                  first = indices[j] == 0;
  4470.  
  4471.               if (first)
  4472.                 {
  4473.                   bounds[whichcomponent * 2] = (real) fval;
  4474.                   bounds[whichcomponent * 2 + 1] = (real) fval;
  4475.                 }
  4476.               else
  4477.                 {
  4478.                   if (fval < bounds[whichcomponent * 2])
  4479.                     bounds[whichcomponent * 2] = (real) fval;
  4480.                   else if (fval > bounds[whichcomponent * 2 + 1])
  4481.                     bounds[whichcomponent * 2 + 1] = (real) fval;
  4482.                 }
  4483.             }
  4484.           else /* read iblank values */
  4485.             {
  4486.               if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4487.                 {
  4488.               if (feof(inFile))
  4489.             {
  4490.                       fprintf(stderr, "Short of data in %s\n", gridptr->name);
  4491.             }
  4492.                   else
  4493.             {
  4494.                       fprintf(stderr, "Error reading %s\n", gridptr->name);
  4495.             }
  4496.                   Free(indices);
  4497.                   Free(bounds);
  4498.                   Free(*formptr);
  4499.                   freedims(ngrid, *dimsptr);
  4500.                   if (gridptr->blank == P3D_BLANK)
  4501.                     {
  4502.                       freeiblank(n + 1, *ipp);
  4503.                     }
  4504.                   fclose(inFile);
  4505.                   return 0;
  4506.                 }
  4507.               putIblank(ival, *ipp, indices, (*dimsptr)[n], ndim, n);
  4508.             }
  4509.  
  4510.           /* adjust the indices and whichcomponent */
  4511.           if (gridptr->permutation == P3D_WHOLE)
  4512.             {
  4513.               for (j = 0; j < ndim; j++)
  4514.                  {
  4515.                    indices[j]++;
  4516.                    if (indices[j] == (*dimsptr)[n][j])
  4517.                      {
  4518.                        indices[j] = 0;
  4519.                        continue;
  4520.                      }
  4521.                    break;
  4522.                  }
  4523.               if (j == ndim) /* all indices are set to zero now */ 
  4524.                 {
  4525.                   whichcomponent++;
  4526.                 }
  4527.             }
  4528.           else  /* P3D_PLANE */
  4529.             {
  4530.               for (j = 0; j < ndim - 1; j++)
  4531.                  {
  4532.                    indices[j]++;
  4533.                    if (indices[j] == (*dimsptr)[n][j])
  4534.                      {
  4535.                        indices[j] = 0;
  4536.                        continue;
  4537.                      }
  4538.                    break;
  4539.                  }
  4540.               if (j == ndim - 1)  
  4541.                 {
  4542.                   whichcomponent++;
  4543.                   if ((gridptr->blank == P3D_BLANK) ? 
  4544.                      (whichcomponent == ndim + 1) : (whichcomponent == ndim))
  4545.                   /* leave a space for reading iblanks */
  4546.                     {
  4547.                       whichcomponent = 0;
  4548.                       indices[j]++;
  4549.                     }
  4550.                 }
  4551.             }
  4552.         } /* for read grid values */
  4553.  
  4554.      (*formptr)[n] = NewCurvilinearDataForm(datasetname, ndim, (*dimsptr)[n],
  4555.                     bounds, dataset);
  4556.  
  4557.     } /* outer for */
  4558.  
  4559.   /* indices and bounds are no longer of any use */
  4560.   Free(indices);
  4561.   Free(bounds);
  4562.   fclose(inFile);
  4563.   return n;
  4564. } /* RdBinXYZ */
  4565.  
  4566. #ifdef PROTO
  4567. static int RdBinSol(int ndim, int isMGrid, int ngrid, int blank,
  4568.                  long **dimsptr, p3d_q *qptr, ObjPtr *formptr, int **ipp)
  4569. #else
  4570. static int RdBinSol(ndim, isMGrid, ngrid, blank, dimsptr, qptr, formptr, ipp)
  4571. int ndim;
  4572. int isMGrid;
  4573. int ngrid;
  4574. int blank;
  4575. long **dimsptr;
  4576. p3d_q *qptr;
  4577. ObjPtr *formptr;
  4578. int **ipp;
  4579. #endif
  4580. {
  4581.   FILE *inFile;
  4582.   int i, j, n;
  4583.   char *name;
  4584.   long *indices;
  4585.   int numsets; /* number of datasets read so far, for multiple datasets in a file*/
  4586.   int namelen; /* the length of dataset name going to be used */
  4587.   int ival;
  4588.   float fval;
  4589.  
  4590.   if (!(inFile = fopen(qptr->name, "r")))
  4591.     {
  4592.       fprintf(stderr, "Unable to open %s, aborted\n", qptr->name);
  4593.       return 0;
  4594.     }
  4595.  
  4596.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  4597.     {
  4598.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  4599.       fclose(inFile);
  4600.       return 0;
  4601.     }
  4602.  
  4603.   /* determine the dataset name from the file name */
  4604.   for (i = (int) strlen(qptr->name) - 1; i >= 0 && qptr->name[i] != '/'; i--)
  4605.      ;
  4606.   if (i == -1 || qptr->name[i] == '/')
  4607.     i++;
  4608.   name = qptr->name + i;
  4609.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  4610.      ;
  4611.  
  4612.   numsets = 0;
  4613.  
  4614.   /* We have to consider multiple datasets here. */
  4615.   do
  4616.     {
  4617.       if (isMGrid)
  4618.     {
  4619.           if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4620.             {
  4621.               if (feof(inFile))
  4622.                 {
  4623.                   if (numsets == 0)
  4624.                     {
  4625.                       fprintf(stderr, "Error in %s: empty file\n", qptr->name);
  4626.                     } /* else it reads to the end of a Mdata file already */
  4627.                   fclose(inFile);
  4628.                   Free(indices);
  4629.                   return (numsets ? ngrid : 0);
  4630.                 }
  4631.               else 
  4632.                 {
  4633.                   fprintf(stderr, "Error reading %s.\n", qptr->name);
  4634.                   fclose(inFile);
  4635.                   Free(indices);
  4636.                   return 0;
  4637.                 }
  4638.             }
  4639.           else if (ival <= 0)
  4640.         {
  4641.               fprintf(stderr, "Error in %s: Invalid NGRID spec.\n", qptr->name);
  4642.               fclose(inFile);
  4643.               Free(indices);
  4644.               return 0;
  4645.         }
  4646.           else
  4647.         {
  4648.           n = ival;
  4649.         }
  4650.         }
  4651.       else
  4652.         {
  4653.           n = 1;
  4654.         }
  4655.  
  4656.       /* if ngrid != this file's ngrid then error */
  4657.       if (ngrid != n)
  4658.         {
  4659.           fprintf(stderr, "Error in %s: NGRID value does not match its grid\n",
  4660.                                                              qptr->name);
  4661.           fclose(inFile);
  4662.           Free(indices);
  4663.           return 0;
  4664.         }
  4665.  
  4666.       /* reading dimensions */
  4667.       for (n = 0; n < ngrid; n++)
  4668.         {
  4669.           for (i = 0; i < ndim; i++)
  4670.             {
  4671.               if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4672.                 {
  4673.                   /* if this is not multiple grid, then the first dimension value 
  4674.              is used to indicate whether there is more datasets or not */
  4675.                   if (!isMGrid && i == 0 && feof(inFile))
  4676.                     {
  4677.                       if (numsets == 0)
  4678.                         {
  4679.                           fprintf(stderr, "Error in %s: empty file\n", qptr->name);
  4680.                         }
  4681.                       fclose(inFile);
  4682.                       Free(indices);
  4683.                       return (numsets ? ngrid : 0);
  4684.                     }
  4685.                   fprintf(stderr, "Error reading dimensions in %s\n", qptr->name);
  4686.                   fclose(inFile);
  4687.                   Free(indices);
  4688.                   return 0;
  4689.                 }
  4690.               /* if dims does not match the dims read in XYZ file then error */
  4691.               if ((long) ival != dimsptr[n][i])
  4692.                 {
  4693.                   fprintf(stderr,
  4694.                   "Dimensions mismatch between %s and its XYZ file\n", qptr->name);
  4695.                   fclose(inFile);
  4696.                   Free(indices);
  4697.                   return 0;
  4698.                 }
  4699.             }
  4700.         }
  4701.  
  4702.       for (n = 0; n < ngrid; n++)
  4703.         {
  4704.           long size, count; /* number of data values promised */
  4705.           ObjPtr density, pressure, dataset;
  4706.           char datasetname[MAXNAMELEN];
  4707.           char densityname[MAXNAMELEN];
  4708.           char pressurename[MAXNAMELEN];
  4709.           double time;
  4710.           int whichcomponent;
  4711.           ObjPtr whichdataset;
  4712.     
  4713.           /* process the FAMACH..TIME line */
  4714.           for (i = 1; i <= 4; i++)
  4715.             {
  4716.               if (1 != fread((void *) &fval, sizeof(float), 1, inFile))
  4717.                 {
  4718.                   fprintf(stderr,
  4719.                      "Error in %s: invalid FSMACH,ALPHA,RE,TIME line\n", qptr->name);
  4720.                   Free(indices);
  4721.                   fclose(inFile);
  4722.                   return 0;
  4723.                 }
  4724.               if (i == 4)
  4725.                 {
  4726.                   time = (double) fval;
  4727.                 }
  4728.             }
  4729.     
  4730.           /* setup datasets */
  4731.           strncpy(datasetname, name, namelen);
  4732.           datasetname[namelen] = '\0';
  4733.           if (ngrid > 1)
  4734.             {
  4735.           /* the ordering starts from 1 */
  4736.               sprintf(datasetname + namelen, "%d", n + 1);
  4737.             }
  4738.           sprintf(densityname, "%s%s", datasetname, " Density");
  4739.           sprintf(pressurename, "%s%s", datasetname, " Pressure");
  4740.           density = NewStructuredDataset(densityname, ndim, dimsptr[n], 0);
  4741.           pressure = NewStructuredDataset(pressurename, ndim, dimsptr[n], 0);
  4742.           dataset = NewStructuredDataset(datasetname, ndim, dimsptr[n], ndim);
  4743.  
  4744.           size = 1;
  4745.           for (i = 0; i < ndim; i++)
  4746.              size *= dimsptr[n][i];
  4747.           size *= ndim + 2;  /* includes density and pressure */
  4748.           for (i = 0; i < ndim; i++)
  4749.              {
  4750.                indices[i] = 0;
  4751.              }
  4752.           whichcomponent = 0;
  4753.           whichdataset = density;
  4754.     
  4755.           /* read in all the q values */
  4756.           for (count = 0; count < size; count++)
  4757.             {
  4758.               int iblank = 1;
  4759.     
  4760.               if (1 != fread((void *) &fval, sizeof(float), 1, inFile))
  4761.                 {
  4762.           if (feof(inFile))
  4763.             {
  4764.                       fprintf(stderr, "Error reading %s, short of data\n", 
  4765.                   qptr->name);
  4766.             }
  4767.                   else
  4768.             {
  4769.                       fprintf(stderr, "Error reading %s\n", qptr->name);
  4770.             }
  4771.                   Free(indices);
  4772.                   fclose(inFile);
  4773.                   return 0;
  4774.                 }
  4775.               /* do check of density and pressure, and iblanks here */
  4776.               if (qptr->check == P3D_CHECK && (whichdataset == pressure ||
  4777.                   whichdataset == density))
  4778.                 {
  4779.                   if (fval < 0.0)
  4780.             {
  4781.               fprintf(stderr, "Negative %s value %f found in Q file %s and is converted to %f\n", whichdataset == pressure ? "pressure" : "density", fval, qptr->name, 0.0);
  4782.                       fval = 0.0;
  4783.                     }
  4784.                 }
  4785.               if (blank == P3D_BLANK)
  4786.                 {
  4787.                    getIblank(&iblank, ipp, indices, dimsptr[n], ndim, n);
  4788.                 }
  4789.     
  4790.               if (!SetCurField(FIELD1, whichdataset))
  4791.                 {
  4792.                   fprintf(stderr, "Unable to SetCurField\n");
  4793.                   Free(indices);
  4794.                   fclose(inFile);
  4795.                   return 0;
  4796.                 }
  4797.     
  4798.               PutFieldComponent(FIELD1, whichcomponent, indices, 
  4799.                           iblank == 0 ? (real) missingData : (real) fval);  
  4800.     
  4801.               /* adjust the indices, whichdataset, and whichcomponent */
  4802.               if (qptr->permutation == P3D_WHOLE)
  4803.                 {
  4804.                   for (j = 0; j < ndim; j++)
  4805.                      {
  4806.                        indices[j]++;
  4807.                        if (indices[j] == dimsptr[n][j])
  4808.                          {
  4809.                            indices[j] = 0;
  4810.                            continue;
  4811.                          }
  4812.                        break;
  4813.                      }
  4814.                   if (j == ndim) /* all indices are set to zero now */ 
  4815.                     {
  4816.                       if (whichdataset == dataset)
  4817.                         {
  4818.                           whichcomponent++;
  4819.                         }
  4820.                       else
  4821.                         {
  4822.                           if (whichdataset == density)
  4823.                             whichdataset = pressure;
  4824.                           else
  4825.                             whichdataset = dataset;
  4826.                         }
  4827.                     }
  4828.                 }
  4829.               else  /* P3D_PLANE */
  4830.                 {
  4831.                   for (j = 0; j < ndim - 1; j++)
  4832.                      {
  4833.                        indices[j]++;
  4834.                        if (indices[j] == dimsptr[n][j])
  4835.                          {
  4836.                            indices[j] = 0;
  4837.                            continue;
  4838.                          }
  4839.                        break;
  4840.                      }
  4841.                   if (j == ndim - 1)  
  4842.                     {
  4843.                       if (whichdataset == density)
  4844.                         whichdataset = pressure;
  4845.                       else if (whichdataset == pressure)
  4846.                         whichdataset = dataset;
  4847.                       else /* whichdataset == dataset */
  4848.                         {
  4849.                           whichcomponent++;
  4850.                           if (whichcomponent == ndim)
  4851.                             {
  4852.                               whichdataset = density;
  4853.                               whichcomponent = 0;
  4854.                               indices[j]++;
  4855.                             }
  4856.                         }
  4857.                     }
  4858.                 }
  4859.             } 
  4860.     
  4861.           SetDatasetForm(density, formptr[n]);
  4862.           SetDatasetForm(pressure, formptr[n]);
  4863.           SetDatasetForm(dataset, formptr[n]);
  4864.           if (qptr->isMDataset)
  4865.             {
  4866.               RegisterTimeDataset(density, (real) time);
  4867.               RegisterTimeDataset(pressure, (real) time);
  4868.               RegisterTimeDataset(dataset, (real) time);
  4869.             }
  4870.           else
  4871.             {
  4872.               RegisterDataset(density);
  4873.               RegisterDataset(pressure);
  4874.               RegisterDataset(dataset);
  4875.             }
  4876.         } /* outer for */
  4877.       numsets++;
  4878.     } while(qptr->isMDataset);
  4879.  
  4880.   Free(indices);
  4881.   fclose(inFile);
  4882.   return ngrid;
  4883. } /* RdBinSol */
  4884.  
  4885. #ifdef PROTO
  4886. static int RdBinFun(int ndim, int isMGrid, int ngrid, int blank,
  4887.                 long **dimsptr, p3d_func *funcptr, ObjPtr *formptr, int **ipp)
  4888. #else
  4889. static int RdBinFun(ndim, isMGrid, ngrid, blank, 
  4890.                               dimsptr, funcptr, formptr, ipp)
  4891. int ndim;
  4892. int isMGrid;
  4893. int ngrid;
  4894. int blank;
  4895. long **dimsptr;
  4896. p3d_func *funcptr;
  4897. ObjPtr *formptr;
  4898. int **ipp;
  4899. #endif
  4900. {
  4901.   FILE *inFile;/* the current function file being processed */
  4902.   int i, j, n;
  4903.   char *name; /* points to the portion of function file name after the last / */
  4904.   long *indices;
  4905.   int *nvar;  /* number of variables of this function */
  4906.   int namelen;
  4907.   int ival;
  4908.   float fval;
  4909.  
  4910.   if (!(inFile = fopen(funcptr->name, "r")))
  4911.     {
  4912.       fprintf(stderr, "Unable to open %s, aborted\n", funcptr->name);
  4913.       return 0;
  4914.     }
  4915.   if (isMGrid)
  4916.     {
  4917.       if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4918.         {
  4919.           fprintf(stderr, "Error reading NGRID value in %s\n", funcptr->name);
  4920.           fclose(inFile);
  4921.           return 0;
  4922.         } 
  4923.       else if (ival <= 0)
  4924.         {
  4925.           fprintf(stderr, "%s format error, invalid NGRID value\n", funcptr->name);
  4926.           fclose(inFile);
  4927.           return 0;
  4928.         } 
  4929.       else
  4930.     {
  4931.       n = ival;
  4932.     }
  4933.     }
  4934.   else
  4935.     {
  4936.       n = 1;
  4937.     }
  4938.   /* if ngrid != this file's ngrid then error */
  4939.   if (ngrid != n)
  4940.     {
  4941.       fprintf(stderr, "Error in %s: NGRID value does not match its grid\n",
  4942.                                                                  funcptr->name);
  4943.       fclose(inFile);
  4944.       return 0;
  4945.     }
  4946.   if (!(indices = (long *) Alloc(ndim * sizeof(long))))
  4947.     {
  4948.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  4949.       fclose(inFile);
  4950.       return 0;
  4951.     }
  4952.   if (!(nvar = (int *) Alloc(ngrid * sizeof(int))))
  4953.     {
  4954.       fprintf(stderr, "Unable to alloc memory, aborted\n");
  4955.       Free(indices);
  4956.       fclose(inFile);
  4957.       return 0;
  4958.     }
  4959.  
  4960.   /* if dims does not match the dims read in XYZ file then error */
  4961.   for (n = 0; n < ngrid; n++)
  4962.     {
  4963.       for (i = 0; i < ndim + 1; i++) /* nvar follows dimension spec's */
  4964.         {
  4965.           if (1 != fread((void *) &ival, sizeof(int), 1, inFile))
  4966.             {
  4967.               fprintf(stderr, "Error reading dim/nvar value in %s\n", funcptr->name);
  4968.               Free(indices);
  4969.               Free(nvar);
  4970.               fclose(inFile);
  4971.               return 0;
  4972.             }
  4973.           if (i < ndim)  /* the value is dimension */
  4974.             {
  4975.               if ((long) ival != dimsptr[n][i])
  4976.                 {
  4977.                   fprintf(stderr, "Dimensions mismatch between %s"
  4978.                                       " and its XYZ file\n", funcptr->name);
  4979.                   Free(indices);
  4980.                   Free(nvar);
  4981.                   fclose(inFile);
  4982.                   return 0;
  4983.                 }
  4984.             }
  4985.           else /* i = ndim, the value is of nvar */
  4986.             {
  4987.               nvar[n] = ival;
  4988.             }
  4989.         }
  4990.     }
  4991.  
  4992.   /* determine the dataset name */
  4993.   for (i = (int) strlen(funcptr->name) - 1; i >= 0 && funcptr->name[i] != '/'; i--)
  4994.      ;
  4995.   if (i == -1 || funcptr->name[i] == '/')
  4996.     i++;
  4997.   name = funcptr->name + i;
  4998.   for (namelen = 0; name[namelen] != '\0' && name[namelen] != '.'; namelen++)
  4999.      ;
  5000.  
  5001.   for (n = 0; n < ngrid; n++)
  5002.     {
  5003.       long size, count; /* number of data values promised */
  5004.       ObjPtr dataset;
  5005.       char datasetname[MAXNAMELEN];
  5006.       int whichcomponent;
  5007.  
  5008.       strncpy(datasetname, name, namelen);
  5009.       datasetname[namelen] = '\0';
  5010.       /* number the datasets only if its number of grids is than 1 */
  5011.       if (ngrid > 1)
  5012.         {
  5013.       /* ordering starts from 1 */
  5014.           sprintf(datasetname + namelen, "%d", n + 1);
  5015.         }
  5016.  
  5017.       size = 1;
  5018.       for (i = 0; i < ndim; i++)
  5019.          size *= dimsptr[n][i];
  5020.       size *= nvar[n];
  5021.       for (i = 0; i < ndim; i++)
  5022.          {
  5023.            indices[i] = 0;
  5024.          }
  5025.  
  5026.       /* setup datasets */
  5027.       dataset = NewStructuredDataset(datasetname, ndim, dimsptr[n],
  5028.                                               nvar[n] == 1 ? 0 : nvar[n]);
  5029.       if (!SetCurField(FIELD1, dataset))
  5030.         {
  5031.           fprintf(stderr, "Unable to SetCurField\n");
  5032.           Free(indices);
  5033.           Free(nvar);
  5034.           fclose(inFile);
  5035.           return 0;
  5036.         }
  5037.       whichcomponent = 0;
  5038.  
  5039.       /* read in all the func values */
  5040.       for (count = 0; count < size; count++)
  5041.         {
  5042.           int iblank = 1;
  5043.      
  5044.           if (1 != fread((void *) &fval, sizeof(float), 1, inFile))
  5045.             {
  5046.           if (feof(inFile))
  5047.         {
  5048.                   fprintf(stderr, "Error reading %s, short of data\n", 
  5049.               funcptr->name);
  5050.         }
  5051.               else
  5052.         {
  5053.                   fprintf(stderr, "Error reading %s\n", funcptr->name);
  5054.                 }
  5055.               Free(indices);
  5056.               Free(nvar);
  5057.               fclose(inFile);
  5058.               return 0;
  5059.             }
  5060.           if (blank == P3D_BLANK)
  5061.             {
  5062.                getIblank(&iblank, ipp, indices, dimsptr[n], ndim, n);
  5063.             }
  5064.     
  5065.           PutFieldComponent(FIELD1, whichcomponent, indices, 
  5066.                   iblank == 0 ? (real) missingData : (real) fval);  
  5067.  
  5068.           /* adjust the indices and whichcomponent */
  5069.           if (funcptr->permutation == P3D_WHOLE)
  5070.             {
  5071.               for (j = 0; j < ndim; j++)
  5072.                 {
  5073.                   indices[j]++;
  5074.                   if (indices[j] == dimsptr[n][j])
  5075.                     {
  5076.                       indices[j] = 0;
  5077.                       continue;
  5078.                     }
  5079.                   break;
  5080.                 }
  5081.               if (j == ndim) /* all indices are set to zero now */
  5082.                 {
  5083.                   whichcomponent++;
  5084.                 }
  5085.             }
  5086.           else  /* P3D_PLANE */
  5087.             {
  5088.               for (j = 0; j < ndim - 1; j++)
  5089.                 {
  5090.                   indices[j]++;
  5091.                   if (indices[j] == dimsptr[n][j])
  5092.                     {
  5093.                       indices[j] = 0;
  5094.                       continue;
  5095.                     }
  5096.                   break;
  5097.                 }
  5098.               if (j == ndim - 1 && ++whichcomponent == nvar[n])
  5099.                 {
  5100.                   indices[j]++;
  5101.                   whichcomponent = 0;
  5102.                 }
  5103.             }
  5104.         }     /* for size */
  5105.  
  5106.       SetDatasetForm(dataset, formptr[n]);
  5107.       RegisterDataset(dataset);
  5108.     } /* outer for */
  5109.  
  5110.   Free(indices);
  5111.   Free(nvar);
  5112.   fclose(inFile);
  5113.   return ngrid;
  5114. } /* RdBinFun */
  5115.  
  5116. /*********************binary p3d reading routines, end *********************/
  5117.  
  5118. #ifdef PROTO
  5119. static void parseDataFiles(int nSets, p3d_set *sets)
  5120. #else
  5121. static void parseDataFiles(nSets, sets)
  5122. int nSets;
  5123. p3d_set *sets;
  5124. #endif
  5125. {
  5126.   int i, j;
  5127.  
  5128. /* If no fortran compiler, then unable to handle unformatted p3d data files. */
  5129.   for (i = 0; i < nSets; i++)
  5130.      {
  5131.        p3d_set *curset = sets + i;
  5132.        ObjPtr *formptr;           /* an array of ngrid dataforms */
  5133.        int **ipp;                 /* points to ngrid arrays of iblank ints */ 
  5134.        long **dimsptr;            /* points to ngrid arrays of dims */
  5135.        int ngrid;
  5136.  
  5137.        if (curset->grid->fileType == P3D_FMT)
  5138.      {
  5139.            ngrid = RdFmtXYZ(curset->ndim, curset->isMGrid, &dimsptr, curset->grid, 
  5140.                 &formptr, &ipp);
  5141.      }
  5142.        else if (curset->grid->fileType == P3D_BIN)
  5143.      {
  5144.            ngrid = RdBinXYZ(curset->ndim, curset->isMGrid, &dimsptr, curset->grid, 
  5145.                 &formptr, &ipp);
  5146.      }
  5147.        else
  5148.      {
  5149. #ifdef FORTRAN
  5150.            ngrid = RdUnfmtXYZ(curset->ndim, curset->isMGrid, &dimsptr, curset->grid, 
  5151.                 &formptr, &ipp);
  5152. #else
  5153.            fprintf(stderr, "Error in %s, %s data reader is not available\n",
  5154.            curset->grid->name, "Unformatted");
  5155.            return;
  5156. #endif
  5157.          }
  5158.        if (!ngrid)
  5159.          {
  5160.            FileFormatError("parseDataFiles", curset->grid->name);
  5161.            return;  /* something wrong with the grid file */
  5162.          }
  5163.  
  5164.        /* read all the q files of this set and construct datasets for them */
  5165.        for (j = 0; j < curset->nQ; j++)
  5166.           {
  5167.         int status;
  5168.  
  5169.             if ((curset->qlist + j)->fileType == P3D_FMT)
  5170.           {
  5171.                 status = RdFmtSol(curset->ndim, curset->isMGrid, ngrid,
  5172.                   curset->grid->blank, dimsptr, curset->qlist + j, formptr, ipp);
  5173.           }
  5174.             else if ((curset->qlist + j)->fileType == P3D_BIN)
  5175.           {
  5176.                 status = RdBinSol(curset->ndim, curset->isMGrid, ngrid,
  5177.                   curset->grid->blank, dimsptr, curset->qlist + j, formptr, ipp);
  5178.           }
  5179.             else
  5180.           {
  5181. #ifdef FORTRAN
  5182.                 status = RdUnfmtSol(curset->ndim, curset->isMGrid, ngrid,
  5183.                       curset->grid->blank, dimsptr, curset->qlist + j, formptr, ipp);
  5184. #else
  5185.                 fprintf(stderr, "Error in %s, %s data reader is not available\n",
  5186.         (curset->qlist + j)->name, "Unformatted");
  5187.                 return;
  5188. #endif
  5189.               }
  5190.             /* return a value of 0 in status if failed */
  5191.             if (!status)
  5192.               {
  5193.                 FileFormatError("Q file", curset->qlist[j].name);
  5194.                 Free(formptr);
  5195.                 freedims(ngrid, dimsptr);
  5196.                 if (curset->grid->blank == P3D_BLANK)
  5197.                   freeiblank(ngrid, ipp);
  5198.                 return;
  5199.               }
  5200.           }
  5201.  
  5202.        for (j = 0; j < curset->nFunc; j++)
  5203.           {
  5204.         int status;
  5205.  
  5206.             if ((curset->funclist + j)->fileType == P3D_FMT)
  5207.           {
  5208.                 status = RdFmtFun(curset->ndim, curset->isMGrid, ngrid,
  5209.                   curset->grid->blank, dimsptr, curset->funclist + j, formptr, ipp);
  5210.           }
  5211.             else if ((curset->funclist + j)->fileType == P3D_BIN)
  5212.           {
  5213.                 status = RdBinFun(curset->ndim, curset->isMGrid, ngrid,
  5214.                   curset->grid->blank, dimsptr, curset->funclist + j, formptr, ipp);
  5215.           }
  5216.             else
  5217.           {
  5218. #ifdef FORTRAN
  5219.                 status = RdUnfmtFun(curset->ndim, curset->isMGrid, ngrid,
  5220.                   curset->grid->blank, dimsptr, curset->funclist + j, formptr, ipp);
  5221. #else
  5222.                 fprintf(stderr, "Error in %s, %s data reader is not available\n",
  5223.         (curset->funclist + j)->name, "Unformatted");
  5224.                 return;
  5225. #endif
  5226.               }
  5227.             /* return a value of 0 in status if failed */
  5228.             if (!status)
  5229.               {
  5230.                 FileFormatError("Function file", curset->funclist[j].name);
  5231.                 Free(formptr);
  5232.                 freedims(ngrid, dimsptr);
  5233.                 if (curset->grid->blank == P3D_BLANK)
  5234.                   freeiblank(ngrid, ipp);
  5235.                 return;
  5236.               }
  5237.           }
  5238.  
  5239.        /* clean up */
  5240.        if (curset->grid->blank == P3D_BLANK)
  5241.          {
  5242.            freeiblank(ngrid, ipp);
  5243.            ipp = NULL;
  5244.          }
  5245.        freedims(ngrid, dimsptr);
  5246.        Free(formptr);
  5247.        dimsptr = NULL;
  5248.        formptr = NULL;
  5249.     } /* outer for */
  5250.   return;
  5251. } /* parseDataFiles */
  5252.  
  5253. #ifdef PROTO
  5254. static ObjPtr ReadP3DFile(char *p3dname)
  5255. #else
  5256. static ObjPtr ReadP3DFile(p3dname)
  5257. char *p3dname;
  5258. #endif
  5259. {
  5260.   int nSets = 0; /* number of p3d sets */
  5261.   p3d_set *sets = NULL; /* pointer to an array of p3d sets */ 
  5262.  
  5263.   /* parse p3d meta file */
  5264.   if (!parseP3DMeta(p3dname, &nSets, &sets))
  5265.     {
  5266.       FileFormatError("ReadP3DFile", p3dname);
  5267.     }
  5268.   else
  5269.     {
  5270.       /* find out if any of the sets using the same XYZ file */
  5271.       compress(&nSets, &sets);
  5272.       parseDataFiles(nSets, sets);
  5273.     }
  5274.  
  5275.   cleanup(nSets, sets);
  5276.   return NULLOBJ;
  5277. } /* ReadP3DFile */
  5278.  
  5279. /*******************************************************************************/
  5280. /*******************************************************************************/
  5281. /*******************************************************************************/
  5282. #ifdef PROTO
  5283. void InitHwuFiles(void)
  5284. #else
  5285. void InitHwuFiles()
  5286. #endif
  5287. {
  5288.   DefineFormat("STF", "stf", ReadSTFFile, 0);
  5289.   DefineFormat("P3D", "p3d", ReadP3DFile, 0);
  5290. }
  5291.  
  5292. #ifdef PROTO
  5293. void KillHwuFiles(void)
  5294. #else
  5295. void KillHwuFiles()
  5296. #endif
  5297. {
  5298. }
  5299.